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1  Introduction 

The  purpose  of  this  report  is  to  assess  the  mathematical  models  currently  avail¬ 
able  for  predicting  the  acoustic  behavior  of  a  seabed  at  frequencies  up  to  about 
500  kHz.  For  unconsolidated  seabeds  there  are  currently  two  rival  models,  the 
general  model  of  a  poroelastic  medium  due  to  Biot  [4]  and  particularized  by 
Stoll  [37]  to  seabeds,  and  the  model  of  Buckingham  [11],  [12],  [13].  The  Biot- 
Stoll  model  treats  a  poroelastic  medium  as  an  elastic  frame  with  interstitial  pore 
fluid.  It  has  two  loss  mechanisms,  intergranular  friction  and  viscous  friction  due 
to  the  motion  of  the  frame  relative  to  the  fluid.  It  depends  upon  thirteen  para¬ 
meters,  which  are  listed  in  Table  1.  The  Buckingham  model  is  derived  from  the 
linearized  Navier-Stokes  equations.  Attenuation  is  due  to  grain-to-grain  friction 
and  strain  hardening  of  the  fluid  between  grains.  It  may  be  regarded  as  a  causal 
version  of  the  standard  elastic  model  of  the  seabed  in  the  sense  that,  like  the 
elastic  model,  it  predicts  wave  attenuations  that  are  approximately  linear  in 
frequency,  but  unlike  the  elastic  model,  it  predicts  wave  speeds  that  are  log¬ 
arithmically  dispersive  as  required  by  the  Kramers-Kronig  causality  relations 
(see  Section  5.4.1). 

Assessing  the  models  requires  addressing  several  questions 

•  To  what  extent  are  the  predictions  of  the  model  in  accord  with  experi¬ 
mental  observations?  In  particular  both  models  were  formulated  prior  to 
Sediment  Acoustics  Experiment  1999  (SAX99)  [41],  [33]  which  is  to  date 
the  most  comprehensive  attempt  to  measure  such  observables  as  wave 
speeds  and  attenuations  for  a  sediment  with  known  parameters. 

•  Are  the  underlying  assumptions  used  in  deriving  the  models  appropri¬ 
ate  for  an  unconsolidated  sediment?  For  instance  does  an  aggregation  of 
uncemented  sand  grains  constitute  a  Hookean  elastic  frame  in  the  sense 
assumed  by  the  Biot  model?  To  the  extent  that  this  is  not  the  case  what 
modifications  are  necessary? 

•  Is  the  Biot  model  unnecessarily  complicated?  The  conventional  model 
Biot-Stoll  model  depends  upon  thirteen  input  parameters,  some  of  which 
are  not  easy  to  measure.  Are  the  viscous  losses  predicted  by  the  Biot 
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model  due  to  the  relative  motion  of  the  fluid  and  frame  important  for 
unconsolidated  seabeds?  If  not  might  the  simpler  Buckingham  model 
which  requires  only  five  inputs  suffice? 

Let  us  begin  by  illustrating  the  different  predictions  the  two  models  make 
in  regards  to  wave  speeds  and  attenuations.  Figures  1  and  2  illustrate  the 
differences.  The  measured  values  of  compressional  wave  speed  and  attenuation 
were  taken  from  Williams  et  al.  [44],  Figure  6.  They  are  a  compilation  of 
measurements  made  during  Sediment  Acoustics  Experiment  1999  (SAX99)  [41], 
[33].  The  parameters  used  in  the  Biot-Stoll  model  to  produce  the  prediction 
labelled  "Best  in  range  fit"  in  the  Figures  are  given  in  Table  1,  Column  4,  and  are 
also  from  [44] .  These  were  chosen  to  give  the  best  fit  to  the  data  using  parameter 
values  that  were  within  the  range  estimates1  in  Table  1,  Column  3.  The  authors 
noted  that  better  agreement  with  the  two  measurements  of  compressional  wave 
speed  below  1  kHz  resulted  from  changing  the  parameters  porosity,  tortuosity 
and  permeability2  to  f3  =  41.5,  a  =  1.12,  k  =  5  x  10~n,  all  of  which  lie  outside 
of  the  estimated  ranges.  These  Biot  predictions  are  labelled  "Best  fit".  As 
indicated  in  Figure  1  the  Biot  model  predicts  that  compressional  wave  speed 
will  be  much  more  dispersive  with  respect  to  frequency  than  the  Buckingham 
model  and  is  in  better  accord  with  the  measured  wave  speeds.  On  the  other 
hand  Figure  2  shows  that  the  Biot  model  underestimates  the  compressional 
wave  attenuation  at  frequencies  above  100  kHz  while  the  Buckingham  model 
is  in  good  agreement  at  all  frequencies.  To  date  no  one  has  offered  a  means 
to  correct  the  Buckingham  model’s  failure  to  predict  the  velocity  dispersion 
observed  in  the  SAX99  exercise  and  in  other  experiments.  There  have  recently 
been  successes  in  improving  the  Biot  model’s  predictions,  however,  and  so  the 
report  will  focus  on  this  model. 


2  Derivation  of  the  equations  of  the  Biot-Stoll 
model 

Biot  developed  his  model  for  a  poroelastic  medium  in  a  series  of  papers  [3] ,  [6] , 
[7],  [5],  [4],  [8]  written  over  the  course  of  more  than  twenty  years.  It  is  a  general 
model  which  has  been  applied  to  areas  of  science  and  engineering  as  diverse  as 
seabeds,  soundproofing  and  ultrasound  propagation  through  bone. 

2.1  Constitutive  equations 

Biot’s  notation  changed  over  the  course  of  the  years.  In  deriving  the  consti¬ 
tutive  equation  and  equations  of  motion  we  follow  Biot  [4]  and  Stoll  [37]  since 
the  notation  in  these  articles  is  the  most  commonly  used  in  seabed  acoustics. 
The  medium  is  regarded  as  an  elastic  frame  with  an  interstitial  pore  fluid.  Two 

1  The  authors’  estimates  of  the  parameter  ranges  are  selective;  they  chose  to  disregard  some 
measurements  of  the  parameters  that  they  deemed  unreliable. 

2  Unless  othrwise  indicated,  the  units  for  physical  parameters  are  meters-kilogram-seconds. 
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Figure  1:  SAX99  measurements  of  compressional  wave  speed.  The  Bucking¬ 
ham  model  predictions  were  computed  from  the  best  fit  parameters  in  Table  1 
assuming  cp  =  1771m  /  s  at  38  kHz,  cs  =  129  m  /  s,a8  =  30  dB/m  at  1kHz. 


Symbol 

Parameter 

Range 

Best  fit 

Pf 

Density  of  the  pore  fluid 

[1020, 1023] 

1023  kg / m3 

Pr 

Density  of  sediment  grains 

[2664, 2716] 

2690  kg / m3 

Kh  =  He  K*h 

Real  frame  bulk  modulus 

4.36  x  107  Pa 

ImA- 

Imag  frame  bulk  modulus 

-2.08  x  10B  Pa 

H  =  Re  fi* 

Real  frame  shear  modulus 

[1.7, 4.3]  x  107 

2.92  x  107  Pa 

Im/i* 

Imag  frame  shear  modulus 

-[1.3, 2.6]  x  10Y 

-1.8  x  10B  Pa 

Kf 

Fluid  bulk  modulus 

[2.37,  2.42]  x  10y 

2.395  x  10y  Pa 

K, 

Grain  bulk  modulus 

[3.2, 4.9]  x  10iU 

3.2  x  101U  Pa 

P 

Porosity 

[0.363, 0.394] 

0.385 

V 

Viscosity  of  pore  fluid 

[0.95,1.15]  x  10~3 

1.05  x  10 3kg/m-s 

k 

Permeability 

[2.1, 4.5]  x  10“n 

2.5  x  10-11  ui2 

a 

Tortuosity 

[1.19,1.57] 

1.35 

a 

Pore  size  parameter 

2.65  x  10-5  m 

Table  1:  Parameter  ranges  for  the  Biot-Stoll  model  measured  in  SAX99.  The 
pore  size  parameter  was  caluculated  from  other  parameters.  The  real  and  imag¬ 
inary  parts  of  the  frame  bulk  modulus  were  taken  from  the  literature. 
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Figure  2:  SAX99  measurements  of  compressional  wave  attenuation.  The  Buck¬ 
ingham  model  predictions  were  computed  from  the  best  fit  parameters  in  Table 
1  assuming  cp  =  1771  m/s  at  38kHz,  cs  =  129 m/s,  as  =  30  dB/m  at  1kHz. 
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displacement  vectors  u(x,y,  z,t)  =  [ux(x,y,z,t),uy(x,y,  z,t),uz(x,y,  z,t)]  and 
U(x,y,z,t)  =  [Ux(x,y,  z,t),Uy(x,y,z,t),Uz(x,y,  z,t)]  track  the  macroscopic 
motion  of  the  frame  and  fluid  respectively  while  the  divergences  e  =  V-u  and 
e  =  V-U  give  the  frame  and  fluid  dilatations.  For  an  isotropic  elastic  body  the 
strain  energy  is  a  function 


W  =  W(h,I2,I3) 

of  the  three  elastic  invariants  (see  Love  [29]) 


I,  =  e 


xx  r  '-"yy  \  '-'zz 


ezz  =  e 


1 


I2  —  eyyezz  +  &XX&ZZ  +  exxeyy  ^  (eyZ  d"  &xz  d"  exy) 

t  1  /  2  2  2  \ 

-*3  —  ^xxCyy&zz  H-  ^  {Cyz&xz^-xy  ^xx^-yz  ^yy^xz  ^zz&xy )  * 

where  the  strains  are  related  to  the  displacements  by 


duy 

&xx  =  o  5  &yy  =  ^  5  = 


dux 
dx 

dux  duv 


duz 
dz 

dur,  du~ 


duy  duz 


^xy  —  ^  n  n  '  o  ’  yz  —  n  ^  o  ’ 

a?/  ox  oz  ox  oz  oy 

For  a  poroelastic  medium  a  strain  energy  function  of  the  form 

w  =  w(h,i2,i3,o 


(i) 


is  posited.  The  additional  quantity  is  the  increment  of  fluid  content 

C  =  /3(e-e)  (2) 

where  ft  denotes  the  porosity  of  the  medium.  For  small  displacements  a  strain 
energy  function  can  be  taken  to  be  a  linear  function  of  the  quadratic  terms  only 

W  =  W(h,I2,  0  =  f  e2  -  2yl2  -  Ce C  +  y  C2- 

The  choice  of  symbols  for  the  coefficients  is  made  simply  to  arrive  at  the  form 
of  the  constitutive  equations  given  in  [4] .  The  six  distinct  aggregate  stresses  on 
an  element  of  the  medium  are  denoted  by 


Fra;  :  Fry :  Fr z  :  ^ ~yy  :  ^ ~yz  :  F2  • 

The  constitutive  equations  are  then  found  by  differentiation  txx  =  ,  rxy  = 

. . .  ,pf  =  where  pf  is  the  pore  fluid  pressure,  to  be 

Fra:  —  He  2/i  i&yy  T  ezz)  C( (  (3) 

xyy  =  He  2/r  {exx  d-  &zz)  C£ 

t~zz  =  He  2/r  ( exx  -f-  ^yy)  C([ 

Fry  l^xy :  Frz  p&xz  >  ^ ~yz  P'^yz 

Pf  =  M( (  —  Ce. 
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From  these  equations  the  modulus  p  is  seen  to  be  the  Lame  shear  coefficient.  In 
the  absence  of  the  two  fluid  coefficients  C  and  M  the  first  four  equations  become 
those  of  an  elastic  solid  with  the  modulus  H  becoming  the  Lame  compressional 
coefficient. 

In  order  to  relate  the  moduli  H,  C  and  M  to  standard  quantities  in  elastic 
and  fluid  media,  Stoll  followed  Biot  and  Willis  [5]  in  considering  two  thought 
experiments: 

•  The  jacketed  test  (Figure  3)  in  which  a  saturated  sample  of  the  porous 
medium  is  placed  in  a  flexible,  impervious  jacket  and  a  constant  external 
pressure  p'  is  applied.  Fluid  is  allowed  to  drain  from  the  sack  in  order  to 
keep  the  fluid  pressure  constant,  whence  it  may  be  normalized  to  pf  =0.  In 
this  case  the  frame  will  exert  an  equilibrating  counter  pressure  Kbe  =  —p' 
where  Kb  is  the  frame’s  bulk  modulus. 

•  The  unjacketed  test  (Figure  4)  in  which  a  sample  of  the  medium  is  im¬ 
mersed  in  fluid  and  a  uniform  pressure  p'  is  applied  to  the  surface  of  the 
fluid.  In  this  case  Pf  =  p'  and  the  counteracting  pressures  are  Kre  = 
Kfe  =  —p'  and  Kfr(j  =  p' ,  where  Kr  is  the  bulk  modulus  of  the  frame 
material,  Kf  is  the  bulk  modulus  of  the  fluid,  and  Kf  r  is  another  modulus. 

In  either  test 


T~XX  -  T~yy  -  TZZ  -  P 

T~xy  —  Txz  TyZ  —  0. 

Upon  substituting  the  first  equation  into  the  first  three  equations  of  (3)  and 
then  adding  the  equations  together  we  obtain 

—p'  =  He  —  4/3/ie  —  C(j.  (4) 

For  the  jacketed  test  we  then  have  from  (4) 

Kb  =  H  -  4/3 n  -  C(/e. 

By  setting  Pf  =  0  in  the  last  equation  in  (3)  we  obtain 

C/e  =  C/M 

whence 

Kb  =  H  -  4/3/x  -  C2/M.  (5) 

For  the  unjacketed  test  we  have  from  (4) 

Kr  =  H  —  4/3/i  +  CKr/Kfr.  (6) 

Also  for  the  unjacketed  test  we  have  from  the  last  equation  in  (3) 

1  =  M/Kfr  +  C/Kr  (7) 
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Solving  (5), (6)  and  (7)  for  H,C  and  M  gives 


(8) 


C  = 
M  = 


Kr(Kr  -  Kb) 
D~Kb 
Kl 


D-Kb 


where 


D  =  Kr(l  +  Kr/Kfr).  (9) 

If  the  sample  undergoes  no  change  in  porosity  during  the  unjacketed  test  then 

and  thus  Kfr  is  related  to  the  bulk  modulus  Kf  of  the  pore  fluid  by 


1 

Wr 


(10) 


in  which  case 

D  =  Kr(l  +  /3(Kr/Kf  —  1)).  (11) 

Observe  that  the  first  two  terms  of  (8)i  are  the  compressional  Lame  coefficient 
of  an  elastic  solid.  The  last  term  then  represents  a  correction  for  the  presence 
of  the  pore  space.  It  goes  to  zero  as  (3  — >  0  and  Kr  — +  Kb. 


2.2  Equations  of  motion 

In  [4]  Biot  derived  equations  of  motion  for  a  poroelastic  medium  from  Lagrange’s 
equations 


w  1  XX 

+ 

3t 

w '  xy 

+ 

drxz  _ 

dx 

dy 

dz 

dTXy 

c)t 

U/yy 

dTyz  _ 

dx 

1 

dy 

1 

dz 

drxz 

+ 

dryz 

+ 

drzz  _ 

dx 

dy 

dz 

9pf  = 

dx 

g PJ  = 

dy 

9pf  = 

dz 

8 


d  /  dT\ 
dt.  V  0ux  ) 
d_  /  8T\ 

dt  \diiy) 

d  /  dT\ 
dt  \diiz) 
d_  (  dT_\  dD_ 

dt  \dwx)  +  dwx 

d_  fdT_\  dD_ 

dt  \dwy)  dwv 

d_  (  dT_\  dD_ 

dt  \dwz)  dwz 


(12) 


Figure  4:  Unjacketed  test. 
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with  the  kinetic  energy  function 


(13) 


(  ■  2  ,  .2  .2 
,  ux  +  u,.  +  uT)  - 


^  -  2/9/  (uxwx  +  iiyWy  +  uzwz)  +  m  {w2x  +  w 2y  +  vy^ 


and  the  Rayleigh  dissipation  function 


D  = 


w2x  +  w2y+  U2 


Here  p  =  (1  —  f3)pr  +  /3pf  is  the  aggregate  density  of  the  medium,  pr  is  the 
density  of  the  frame  material,  p/  is  the  density  of  the  pore  fluid,  and  m  and  d 
are  mass  and  dissipation  parameters  which  will  be  discussed  shortly.  The  vector 


w  =  fi(u  —  U) 


measures  the  relative  displacement  of  the  frame  and  fluid.  The  dots  over  the 
components  of  u  and  w  denote  partial  derivatives  with  respect  to  time.  When 
the  expressions  for  T  and  D  are  substituted  into  (12)  the  following  equations  of 
motion  result 

o  2 

pS72u  +  V[(H  -  p)e-CC,\  =  -^(pu-  pfw)  (14) 

V[Ce  -  MQ  =  ( Pfu  -  mw)  ~  d^f  w- 

For  uniform  pore  size  and  a  flow  parallel  to  the  pore  direction,  the  mass 
parameter  m  would  be  pf/0,  so  that  in  (14)2  p/u  —  mw  =  PfU,  but  given  that 
not  all  of  the  fluid  will  flow  in  the  direction  of  the  macroscopic  pressure  gradient 
due  to  the  tortuosity  of  the  pore  space,  Stoll  [37]  takes  the  parameter  m  to  be 

apf 

m  =  ,  oi  >  1 


where  a  is  referred  to  as  the  tortuosity .  Stoll  states  that  for  pores  with  all 
possible  orientations  the  value  of  a  would  be  3,  but  he  gives  no  explanation  or 
citation. 

Finally  note  that  the  dissipation  term  in  (14)2  has  the  dimensions  of  a  pres¬ 
sure  gradient.  According  to  Darcy’s  law  for  steady  laminar  (Poiseuille)  flow 
through  a  porous  medium  the  average  seepage  velocity  (V)  due  to  a  pressure 
gradient  ||  is 


(V)  =  ~ 


k  dp 


/3r]  dx 

where  (3  and  k  are  the  porosity  and  permeability  of  the  medium. 
— /?  (V)  with  jjj  W  in  (14)2,  suggests  that  the  parameter  d  should  be 


(15) 

Identifying 


d  = 


V 

k' 


This  is  the  case  for  only  Poiseuille  flow,  however.  The  modification  required  for 
oscillatory  motion  is  discussed  in  the  next  section. 

3It  is  also  referred  to  as  the  structure  factor  or  constant  or  the  virtual  mass  coefficient. 
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2.3  The  general  form  of  the  dissipation  parameter 

To  find  the  form  of  the  dissipation  parameter  d  in  (14)2  we  need  to  calculate 
the  ratio  of  the  force  per  unit  volume  exerted  on  the  pore  walls  to  the  aver¬ 
age  velocity  of  the  fluid  flowing  through  the  pore  space.  Biot  [7]  solved  the 
problem  of  fluid  motion  in  a  duct  for  two  cases  that  he  regarded  as  opposing 
extremes,  a  cylindrical  duct  and  a  duct  between  two  parallel  plates,  the  latter 
approximately  representing  the  case  of  flat  ovaloid  pores.  Consider  a  flow  in 
the  ^-direction  for  each  case.  For  a  cylinder  of  radius  a  with  spatial  coordinates 
(x,  r,  0)  considerations  of  symmetry  dictate  that  the  velocity  of  the  flow  V  will 
depend  only  on  r  and  time  t  and  the  stress  r  on  the  wall  of  the  cylinder  only 
on  t.  The  dissipation  parameter  d  will  take  the  form 


d  = 


2'KarAx 

7va2Ax 


HV) 


2  r 

(3a  (V) 


(16) 


where  ira2Ax  is  the  volume  of  an  arbitrary  element  of  the  cylindrical  duct  and 
( V )  denotes  the  average  seepage  velocity.  On  the  other  hand  for  two  plates  of 
separation  2ai  with  spatial  coordinates  ( x ,  y,  z ),  — ai  <  y  <  a\  considerations  of 
symmetry  dictate  that  the  stress  r  is  the  same  on  the  top  and  bottom  plates 
whence  for  an  arbitrary  element  of  the  duct 


d  = 


2tAxAz 

AxAyAz 


(3(V) 


T 

(3a i  (V) ' 


(17) 


The  position  of  the  fluid  U(r,t )  is  governed  by  the  linearized  Navier-Stokes 
equation 

d2U  dp  .  fdU N 
Pf~di2 


— s+”v- '  « 


Here  y  is  fluid  viscosity  and  the  pressure  gradient  is  assumed  to  depend  upon 
time,  but  to  be  independent  of  the  spatial  coordinates.  The  position  u(t)  of  the 
walls  of  the  duct  is  also  assumed  to  be  independent  of  x.  The  relative  velocity 
V(r,t)  =  Mi  ( U(r,t )  —  u(t))  then  satisfies 


dV  dp  d2u 


(18) 


Steady,  laminar  flow  in  which  the  inertial  terms  in  (18)  are  neglected 

is  termed  Poiseuille  flow.  For  the  case  of  cylindrical  ducts  the  solution  of  the 
differential  equation 


is 


fd2V  1  dV\  _  dp 


V  dr 2 


r  dr  )  dxV  ^  ° 


T  Tf  \  1  (3p  /  2  2\ 


11 


and  the  average  velocity  is  given  by 


{v)=L^pjLfa  r  (r2  _ a2) rdedr = 

4?7  dx  TTCL2  Jo  Jo  y  '  877 


a2  <9p 
dx' 


(19) 


4t;  da;  ?ra2  70  70 
Comparison  of  (19)  with  Darcy’s  law  (15)  gives 

k  =  ^  (20) 

for  the  permeability  of  a  cylindrical  duct.  The  stress  at  the  wall  of  the  duct  is 

r  =  -nV'(a)  =  -“-fx  (21) 

and  thus  from  (16) 


d  = 


2  r  877  77 

/3a  (C) 


For  the  case  of  Poiseuille  flow  in  a  pore  duct  between  parallel  plates  equation 
(18)  gives 

d2V  dp 

=  ^,C(±a!)  =0 


dy 2 


dx 1 


which  has  solution 


The  average  velocity  is  then 


/T/\  1  3?  1  r  (  2  2\  ,  1  9p  2 


277  5a;  2ai  J_a 
From  Darcy’s  law  (15)  it  follows  that 


jfe  = 


(22) 


The  stress  on  the  wall  is 


and  thus  from  (17) 


r  =  —rjV'{a\)  =  -at 


dp 

dx 


d  = 


377  77 


(3ai  (V)  @a\ 


If  the  duct  is  assumed  to  be  undergoing  time-harmonic  oscillations  V (x,  y,  z,  t) 


V(x,y,z)eiu;t,  |f  (t)  =  |f  e^Mt)  = 
Xelut,  then  (18)  becomes 


m  =  (-i£+^“) 


yV2V  —  iuipjV  =  —pfX. 


(23) 
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Note  that  u>  =  0  corresponds  to  Poiseuille  flow.  For  the  case  of  a  cylindrical 
duct  (23)  has  solution 


V(r) 


X_  (  _  J0  ( i3/2nr/a )  \ 
ioj  y  Jo  (z3/2k)  J 


where 


The  average  velocity  is  then 


ee  X(u)a. 


.  .  2  ra  .  .  X  (  2ViJi  U3/2k)  \ 

(V)  =  “2  /  V{r)rdr=  —  il  + -  ,.3/9  v 

a  Jo  y  kJ0  (z3/“k)  J 


(24) 


X  (  2i  ,  ,  .  ,  s 
-  1  +  -Tk  25 

!W  \  K 


where 


The  stress  at  the  wall  is 


TM  .  4^ 


ViJo(i3/2K) 


t  =  —rjV'{a)  =  —  ??T 


X  z3/2kJi  (z3/2k)  tjkX 


iui  aJ0  (z3/2k)  zwa 


r(«) 


(26) 


and  thus 


where 


d  = 


2r 


Pa{V) 


8?7 

/3a2 


T(K) 


V 

k 


F(k) 


F(n) 


1  kT(k) 

4  1  +  2fT(/s) ' 


Figure  5  shows  the  function  F(k).  As  can  be  seen  F(k)  — >  1  as  k  — >  0,  thus 
giving  d  =  r]/k  for  Poiseuille  flow  in  the  low  frequency  limit. 

For  the  case  of  flow  between  two  plates  of  separation  2ai  we  obtain 


d  =  ^Fi(ki) 


where 


Ki  =  x(w)ai  ,Fi(k) 


1  tanlr 
3  1  +  ^  tanlr  vJk 


Observing  that  the  curves  Re.Fi(/ci)  and  ImFjsi)  were  almost  indistinguish¬ 
able  from  the  curves  ReF(n)  and  IhiFJk)  with  a  =  4a±/3  in  (24),  Biot  conjec¬ 
tured  that  it  might  be  possible  to  compensate  for  the  deviation  from  Poiseuille 
flow  with  increasing  frequency  by  using  the  complex-valued  dynamic  viscosity 
r]F(n)  with  the  pore  size  parameter  a  in  (24)  being  chosen  appropriately  for  the 
dimensions  and  geometry  of  the  pores.  He  did  not  give  much  guidance  as  to 
how  to  measure  the  pore  size  parameter  other  than  to  suggest  is  might  be  ascer¬ 
tainable  from  experimental  data  on  dispersion  of  wave  speed  and  attenuation 
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3 


Figure  5:  Biot’s  function  F(k). 


with  respect  to  frequency.  The  dissipation  parameter  d  in  (14)2,  then  is  given 

by 

d  =  j:F  (xM«) 

for  any  pore  shape  for  which  the  pore  size  parameter  is  ascertainable.  When  the 
time-harmonic  factor  is  e~iult,  rather  than  elwt,  then  d  is  replaced  by  its  complex 
conjugate  d.  Figure  6  shows  the  effect  that  the  presence  of  the  factor  F(k)  has 
on  the  Biot  model’s  prediction  for  compressional  wave  speed  and  attenuation 
when  the  parameters  of  Table  1,  Column  4  are  used. 

2.4  Incorporation  of  viscoelastic  intergranular  effects 

The  only  loss  mechanism  considered  thus  far  is  viscous  loss  due  to  the  relative 
motion  of  the  frame  and  the  pore  fluid.  Clearly  other  loss  mechanisms  such  as 
intergranular  friction  may  be  of  significance.  Following  the  suggestion  of  Stoll 
[37]  (see  also  [39])  most  authors  have  attempted  to  incorporate  such  losses  by 
giving  the  moduli  /i  and  Kb  small  imaginary  parts  which  are  independent  of 
frequency  (see  Section  3.5).  However,  as  noted  by  Turgut  [42],  this  model  is  not 
consistent  with  the  Kramers-Kronig  causality  relations.  Moreover  recent  work 
by  Chotiros  and  Isakson  [17]  suggests  that  a  more  sophisticated  approach  to 
modelling  microscopic  granular  effects  may  rectify  the  Biot  models  poor  predic¬ 
tions  for  compressional  wave  attenuation  at  high  frequencies  that  is  illustrated 
in  Figure  2. 
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Figure  6:  The  blue  line  is  the  Biot  model  prediction  when  d  =  r]F(n)/k.  The 
red  line  is  the  prediction  when  d  =  r)/k.  The  upper  graph  is  compressional  wave 
speed  and  the  lower  graph  is  compressional  wave  attenuation.  The  Biot-Stoll 
parameters  are  the  best  fit  parameters  of  Table  1. 
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FlC.  4.  liiuztrutkn  of  irlual  jrrain  and  fluid  gup  (.’onfifurttion!. 

Figure  7:  From  [8]. 


In  [4]  Biot  gave  the  general  form 


H* 

M* 

C* 


4  *  k*  +  7  -  8 
3^  k*(7  +  S)  —  52 

1 

7  +  S  —  52/n* 

1  -  S/k* 

7  +  1 5  —  52/n*  ' 


Algebraically  these  equations  are  equations  (8)  rendered  in  terms  of  the  com¬ 
pressibilities 

«*  =  K^,5  =  K~\c  =  KJ \  7  =  f3(c  —  S) , 

however  Biot  suggests  that  these  moduli  be  regarded  as  operators  calculated 
from  the  operators  n*  and  p*  which  are  chosen  appropriately  for  the  phenom¬ 
enon  to  be  modelled.  We  will  discuss  three  choices  of  these  operators  which  are 
relevant  to  the  work  of  Chotiros  and  Isakson. 

For  fluid  filled  cracks  and  gaps  (Figure  7)  Biot  [8],  assuming  that  the  normal 
pressure  P  on  the  gap  and  the  change  in  gap  width  Ah  are  related  by  P  =  ZAh, 
arrives  at  the  formula  (Figure  8) 


k  =  k0  + 


1/ki 


(27) 


For  a  fluid-filled  gap  Biot  finds  that  pore  fluid  pressure  satisfies  the  diffusion 
equation 

d2pf  12?7  cdpf 

dx 2  h2  dt 

where  rj  is  viscosity  and  h  is  the  width  of  the  gap  (Figure  9).  Thus,  assuming 
Pf( 0,  t)  =  Pf(D,  t)  =  0,pf(x,  0)  =  const, 


OO 

Pf  =  Bn  exp 

71—1 


12?7c  /  2n  —  1  \ 


sin 


2  n-  1 
D 


-7 TX 
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Fig.  7.  Fluid-aaturatert  crack  At  in  a  solid,  and  it* 
equivalent  block  diagram. 


Figure  8:  From  [8]. 


I®.  5.  Simplified  model  for  a  oinprcssible  fluid  gap 
between  uflid  boundaries.. 


Figure  9:  From  [8] 
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Ftfi.  6.  Maiwell  ebrants 
in  parallel  equivalent  to  the 
flwd  gap  of  Fig,  S 


From  [8]. 


This  suggested  to  Biot  that  a  fluid-filled  gap  should  be  modelled  as  a  series  of 
Maxwell  elements  in  parallel  (Figure  10),  giving 


P  +  rn 


where  p  is  in  general  the  time-derivative  operator  d/dt  and  for  the  specific  case 
of  time-lrarmonic  oscillations  p  =  —  ico  and  the  relaxation  frequencies  are 


rn 


Upc  (2n-l  V 

\  D  / 


Thus  Biot’s  jacketed  compressibility  operator  for  a  fluid  filled  crack  is 


k  =  K0  + 


-I  Y^°°  lOJdr, 

1/^1  Z^n= 1  rn-it 


(28) 


in  the  time-harmonic  case. 

Biot  used  the  following  general  form  for  the  shear  modulus 

r°°  p 

P*=  — ; — p(r)dr  +  p0  +  mp. 

Jo  P  +  r 

For  solid  dissipation  due  to  shear  Biot  [8]  suggests  the  representation 


P*  =  Po  +  a(—icu)s,  0  <  s  <  1 


(29) 


which  follows  from  assuming  the  relaxation  spectrum 

p(r)  =  — rs_1sins7r. 

7 r 

Finally  Chotiros  and  Isakson  [17]  suggest  that  the  presence  of  air  bubbles  in 
the  pore  space  might  reduce  the  effective  fluid  bulk  modulus  (see  Section  3.3) 
and  thereby  explain  the  lower  than  predicted  wave  speeds  below  1  kHz  for  the 
SAX99  sediment  that  are  shown  in  Figure  1  without  assuming  a  higher  than 
measured  porosity.  In  Biot  [8]  posits  that  the  presence  of  air  bubbles  in  the 
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Figure  11:  From  [8] 


pore  fluid  can  be  modelled  by  assuming  the  fluid  compressibility  is  given  by  the 
operator  (Figure  11) 


c*  =  c0  + 


1 

mauj2  —  ibaw 


(30) 


Biot  points  out  that  if  the  damping  constant  ba  is  small  relative  to  the  spring 
constant  ka,  then  there  will  be  a  resonant  frequency. 


2.5  Extensions  and  continuations  of  Biot’s  formulation 

Biot  recognized  that  his  analysis  assumed  that  the  pore  size  distribution  is  nar¬ 
row  and  that  some  generalization  would  be  required  for  a  material  in  which 
pore  sizes  were  more  widely  distributed.  In  this  section  we  present  two  subse¬ 
quent  extensions  of  Biot’s  formulation,  the  formula  for  the  pore  size  parameter 
developed  by  Johnson,  Koplik  and  Dashen  [25]  and  the  generalization  of  the 
viscosity  correction  factor  F  of  Yamamoto  and  Turgut  [45] .  In  this  section  we 
also  describe  the  work  of  Chotiros  and  Isakson  [17]  on  intergranular  viscoelastic 
mechanisms,  which  represents  a  continuation  of  the  work  described  in  Section 
2.4. 


2.5.1  The  pore  size  parameter  for  distributed  pore  sizes 

Johnson,  Koplik  and  Dashen  [25]  argue  that  for  pores  represented  by  non¬ 
intersecting  canted  cylindrical  tubes  of  radius  R 

~Jr2~  = 

where  is  the  limit  of  tortuosity  as  frequency  goes  to  infinity  and  ko  is  the 
limit  of  permeability  as  frequency  goes  to  zero  (in  their  formulation  a  and  k  are 
frequency  dependent)  and  offer  the  heuristic  generalization 


8aoc  &0 

/3A2 


(31) 
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where  «  1  and  A  is  defined  by 


2_fSp\uP\2dA 
A  JVp\up\dV 

where  up  is  the  microscopic  potential  flow  field,  Sp  is  the  pore-grain  surface  and 
Vp  is  the  pore  space  volume.  For  cylindrical  pores  the  flow  field  is  constant  in 
which  case  (32)  gives  A  =  R,  the  pore  radius.  Upon  conducting  simulations  upon 
cubical  lattices  in  which  the  pore  radii  were  distributed  according  to  various 
probability  distributions  they  concluded  that  the  approximation  «  1  "may 
be  violated  by  a  factor  of  2  but  not  apparently  by  a  factor  of  10."  If  we  take 
<^  =  1  and  identify  A  with  the  pore  size  parameter,  and  a  and  k  in  the  Biot-Stoll 
parameter  set  with  a<x,  and  ko  (see  below)  then  we  arrive  at  the  formula 

«  =  (33, 

for  the  pore  size  parameter  for  cylindrical  ducts  of  varying  radii.  For  flat  pores 
the  8  is  replaced  by  3.  While  the  reasoning  is  ad  hoc,  the  formula  has  found 
favor.  It  was  used  to  determine  pore  size  for  seabed  for  SAX99  by  [44]  (see 
Table  1).  It  is  also  used  in  the  application  of  the  Biot  model  to  cancellous  bone 
[20]. 

As  mentioned  above  [25]  assume  that  permeability  and  tortuosity  are  dy¬ 
namic.  An  alternative  formulation  of  Biot’s  equations  due  to  Hovem  and  In¬ 
gram  [24]  gives  a  dynamic  permeability  and  tortuosity.  In  (14)2  take  m  = 
cXooP///?,  d  =  r]F /ko  =  i)(Fr—iFi)/k o  and  assume  time-harmonic  motion  e(x,  y,  z ,  t) 
e(x,  y,  z)e~lut, _  The  result  is 

V[Ce  —  MQ  =  —oj2pfU  + 

P  \  Pfk0u)  k0 

and  thus  the  dynamic  permeability  and  tortuosity  are  given  by 

,/  s  ko  ,  s.  .  vPFi 

fc(w)  =  —  ,a(w)  =  aoc  H - - — .  (34) 

tr  Pfkoto 

Since  F)  is  asymptotically  proportional  to  t/lo,  the  second  term  in  the  formula 
for  a(co)  goes  to  zero.  Likewise  Fr( 0)  =  1  This  justifies  the  identification  of 
and  ko  with  the  Biot-Stoll  tortuosity  parameter  a  and  permeability  parameter 
k  in  (31).  Figure  12  shows  fc(w)  and  a(u j)  for  the  sand  of  SAX99. 

2.5.2  The  viscosity  correction  factor  for  distributed  pore  sizes 

Yamamoto  and  Turgut  [45]  assume  that  the  pore  sizes  are  distributed  according 
to  a  density  function  e(r)  with 

e(r)dr  =  1. 
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Figure  12:  Dynamic  permeability  and  poresize  given  by  formulae  (34)  for  the 
parameters  of  Table  1,  Column  4. 


In  their  scheme  permeability  is  related  to  the  pore  size  distribution  by 


for  cylindrical  and  flat  pores  respectively.  Note  that  this  yields  the  formulas 
(20)  and  (22)  for  the  case  of  a  single  pore  size  e(r)  =  S(r  —  a). 

For  cylindrical  pore  ducts  the  force  on  the  walls  per  unit  volume  is  (see 
(16), (24), (26)) 


rMr.w)  MJ  f2lHITMnMJ 

/  - e(r)dr  =  / - : - 1  (y(w)r)  e(r)ar, 

Jo  r  Jo  r 


while  the  average  velocity  from  (25)  is 


X 


2  i 


( V )  (r,ui)e(r )dr  =  J  —  T(x(u)r)  )  e(r)dr. 


Using  (16)  this  gives 


2^X'(^)  ir  rT  U(u)r)  e(r)dr 
I3  /o°°  (x  +  FuWT(xMr))  e(r)dr 


(36) 
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Observe  that  for  the  case  of  Poiseuille  flow  we  have  from  (21), (19)  and  (35) 

/0°°  2^1e(r)dr  =  877  f0°°  e(r) dr  = 

^/0°°  (P)  (r,  w)e(r)dr  /?  /0°°  r2e(r)dr  kc 

and  thus  the  expression  for  d  has  the  expected  limit  as  cu  — ►  0. 

Yamamoto  and  Turgut  suggest  a  (^-normal  distribution  for  pore  size 

f{4>)  =  ^  exp  ^  2^')  j  ,  -00  <(j)<  00 .  (37) 

Phi  units  are  defined  by  (f>  =  —  log2  r  where  r  is  measured  in  millimeters.  In 
this  case  (36)  becomes 

d_  2?7xM  fZo  VmT  (X(w)r  f  (<j))  dcf) 

PS-o o(1  +  r(fa(M)Tfe^)r  M))  /  M  # 

For  d  =  dr  —  idi  this  gives  the  equations 

a(o>)  =  aoo  +  k(uj)  =  (38) 

ujpf  ar 

for  the  dynamic  tortuosity  and  permeability. 

For  the  </>-normal  distribution  (37)  the  integral  in  (35)  can  be  computed 
exactly,  giving  the  permeability 

k  =  —  2-2^exp(2(CTrf,ln2)2)10-6m2 
np 

where  np  =  3  or  8  for  the  case  of  flat  or  cylindrical  pores  respectively.  Thus  for 
given  k  and 

H  =  ~2lhln(10/3nP  exP (-2(o>ln2)2)^  •  (39) 

Since  the  pore  size  parameter  is  no  longer  used,  there  is  not  an  increase  in  the 
number  of  parameters  from  the  basic  Biot-Stoll  model,  as  long  as  one  accepts 
the  formulas  (35).  However  values  for  the  parameters  and  cr^  are  not  readily 
available.  Yamamoto  and  Turgut  [45]  cite  only  one  instance  of  direct  measure¬ 
ment  of  pore  size  distribution,  "Eagle  River  sand",  for  which  n<j>  ~  5.75,  a $  s=s  1.5 
phi  units. 

Figure  13  shows  how  varying  the  standard  deviation  cr^  of  the  pore  size 
distribution  affects  the  dynamic  permeability  and  tortuosity  computed  from  (38) 
(cf.  Figure  12).  The  other  parameters  used  were  the  those  of  Table  1,  Column  4. 
For  a  fixed  value  of  permeability  increasing  decreases  dynamic  permeability 
and  increases  dynamic  tortuosity  with  the  changes  being  most  pronounced  at 
high  frequencies  for  permeability  and  low  frequencies  for  tortuosity.  Figure 
14  shows  that  the  effect  of  increasing  a<j,  on  compressional  wave  speed  and 
attenuation  is  to  make  the  sediment  behave  more  like  a  viscoelastic  solid,  that 
is,  it  is  less  dispersive  with  respect  to  wave  speed  and  closer  to  log-linear  with 
resect  to  attenuation. 
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Figure  13:  Effect  on  dynamic  permeability  and  tortuosity  of  distributed  pore 
sizes.  The  permeability  was  held  constant,  and  the  mean  pore  size  was  computed 
from  (39). 
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Figure  14:  Effect  on  compressional  wave  speed  and  and  attenuation  of  distrib¬ 
uted  pore  sizes.  The  permeability  was  held  constant,  and  the  mean  pore  size 
was  computed  from  (39). 


24 


2.5.3  Incorporation  of  "squirt  flow 


Biot’s  modelling  of  viscoelastic  forces  at  the  granular  level  is  discussed  above 
in  Section  2.4.  Recently  Chotiros  and  Isakson  [17]  have  proposed  frequency  de¬ 
pendent  forms  for  the  moduli  and  fi*  based  on  similar  considerations.  They 
take  the  region  between  two  sand  grains  to  consist  of  an  area  of  solid  contact 
and  a  gap  into  and  out  of  which  fluid  can  flow  when  the  frame  is  subjected  to 
expansion  and  compression  (Figure  15).  Following  Dvorkin  and  Nur  [18],  they 
term  this  "squirt  flow". The  solid  contact  is  modelled  as  a  spring  of  modulus  Kc 
and  the  effect  of  the  flow  in  the  gap,  which  will  depend  on  the  compressibility 
of  the  fluid  and  the  resistance  to  radial  expansion  normal  to  the  direction  of  the 
force  is  treated  as  a  serial  spring-dashpot  (Maxwell)  element  with  modulus  Ky 
and  damping  constant  r]y.  The  two  mechanisms  are  assumed  to  act  in  parallel 
resulting  in  what  is  referred  to  in  viscoelastic  theory  as  the  standard  linear  el¬ 
ement  (Figure  16(c)).  Under  a  shearing  force  the  absence  of  fluid  loading  for 
shear  causes  the  viscoelastic  element  to  simplify  to  a  Kelvin- Voigt  element  with 
modulus  Gc  and  damping  constant  ijs  (Figure  16(d)).  The  effective  moduli  for 
the  two  types  of  elements  undergoing  time-harmonic  oscillations  e_wt  are  (cf. 
Gittus  [19],  for  instance) 


K 

K*  =  K  A _ v 

b  c+l  +  i(Wfc/o;) 


H*  =  Gc  (  1  -  i  — 

V 

Ky  GC 

Uk  =  — ,  £*v  =  — 
Vv  Vs 


(40) 


As  indicated  in  (40),  Chotiros  and  Isakson  identify  these  effective  moduli  with 
the  Biot-Stoll  parameters  and  n* .  Their  argument  for  this  is  that  the  com¬ 
pressibility  of  the  frame  minerals  is  much  less  than  that  of  the  fluid  and  thus 
storage  and  loss  mechanisms  modelled  by  the  two  viscoelastic  elements  will  pre¬ 
dominate.  They  term  the  resulting  modified  Biot  model  the  "Biot  model  with 
grain  contact  squirt  flow  and  shear  drag"  (BICSQS). 

The  zero-frequency  values  Kc  and  Gc  of  the  two  moduli  are  assumed  to  be 
related  by  an  equation  of  elasticity4 


l  +  v 
1  -  2v 


where  v  is  the  Poisson  ratio  and  thus  the  four  independent  Biot-Stoll  parameters 
Re  K£ ,  Im  ,  Re  /x* ,  Im  /x*  are  replaced  by  the  five  parameters  Gc ,  Ky ,  v  and  the 
two  relaxation  frequencies 

r  _  Wk  j,  _  ^/i 

Jk  -  2^’U  -  2^- 

4  The  of  use  such  elastic  relations  has  been  challenged  by  Chotiros  [15]  on  the  grounds  that 
an  unconsolidated  sediment  is  not  a  Hookean  elastic  solid,  as  is  assumed  in  the  Biot  model. 
In  this  instance  the  parameter  Kc  turns  out  to  have  too  little  effect  on  the  model’s  predictions 
to  test  the  validity  of  this  criticism,  however. 
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FIG.  4.  Grain  contact  models  of  (a)  compression,  (b)  shear,  and  their 
equivalents  (c)  and  (d). 


Figure  16:  From  [17]. 


Comparing  Chotiros  and  Isakson’s  equation  (40) i  with  Biot’s  equation  (28), 
we  see  that  they  are  mathematically  identical  upon  setting  the  compressibility 
Ko  and  the  higher  mode  moduli  an,n  >  1  to  zero  and  making  the  identifications 
Ki  =  K~x,ai  =  Ky ,  n  =  tjJk-  Thus  while  Biot  does  not  state  explicitly  what  the 
compressibilities  Kq  and  /sq  mean  physically,  if  we  assume  that  Kq  represents  the 
compressibility  of  the  sand  grains  and  «q  represents  the  compressibility  of  the 
solid  contact  region  and  assume  as  Chotiros  and  Isakson  do  that  «q  «o  then 
the  two  models  are  both  physically  and  mathematically  similar.  Biot’s  model 
does  not  specifically  mention  resistance  to  radial  flow  as  a  damping  mechanism, 
however.  Likewise  Chotiros  and  Isakson’s  equation  for  the  shear  modulus  (40)2 
is  Biot’s  equation  (29)  with  /io  =  Gc,s  s=s  l,o  =  ^lqoj~s,  though  the  assumed 
physical  mechanisms  are  different. 

Because  the  speed  of  shear  waves  is  controlled  almost  exclusively  by  the  real 
part  of  n*  and  the  parameters  affecting  the  density  of  the  seabed,  the  modulus 
Gc  can  be  computed  from  a  measured  value  of  the  shear  wave  speed,  however 
there  is  no  obvious  way  of  determining  the  remaining  parameters  other  than 
using  them  to  fit  observed  data  such  as  that  displayed  in  Figures  1  and  2.  Figures 
17  and  18  indicate  the  influence  of  the  four  parameters  Ky,v,fk  and  /M  on 
compressional  wave  speed  and  attenuation  in  the  case  of  the  SAX99  sand  (Table 
1,  Column  4)  when  they  are  varied  about  Ky  =  1  GPa,  v  =  0.15,  fk  =  10  kHz 
and  =  50  kHz.  Poisson  ratios  in  the  range  [0.05,0.25]  have  little  effect  on 
either  speed  or  attenuation  in  the  frequency  range  100  Hz  to  1  MHz.  None  of  the 
four  parameters  have  much  effect  below  1  kHz.  The  shear  relaxation  frequency 
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Figure  17:  The  effect  on  compressional  wave  speed  of  varying  four  parameters 
in  the  BICSQS  model:  Ky  =  1  ±  0.2  GPa,  fk  =  10  ±  5  kHz,  ffl  =  50  ±  25  kHz, 
v  =  0.15  ±  0.1.  Green,  blue  and  red  curves  correspond  to  the  lower,  middle  and 
upper  values  respectively. 


ffi  is  influential  only  on  attenuation  in  the  upper  range  of  frequencies.  The 
bulk  modulus  Ky  and  the  bulk  relaxation  frequency  influence  wave  speed  in  the 
middle  to  high  range  of  frequencies  and  attenuation  in  the  middle  range. 

As  discussed  above,  Biot  proposed  models  of  viscoelastic  loss  that  are,  at 
least  mathematically,  generalizations  of  the  BICSQS  model.  Figures  19  and  20 
indicate  the  influence  of  the  other  parameters  and  of  air  bubble  resonance.  It 
was  assumed  in  (28)  that 

/2  I  +  jA"1 
Ko  +  Ki-{3Gc1^2v)  ■ 

Increasing  the  compressibility  kq  lowered  the  compressional  wave  speed  in  the 
high  frequency  range  and  thus  countered  the  effect  of  k±  =  Ky1.  Likewise 
lowering  s  in  equation  (29)  counteracted  the  effect  on  attenuation  of  increasing 
fn,  however  it  also  caused  wave  speed  to  rise  slightly  at  high  frequencies.  Thus 
if  one  is  seeking  to  change  the  predictions  of  the  Biot  model  in  order  to  produce 
better  agreement  with  data,  the  approximations  kq  «  0,  s  ss  1  which  lead  to 
the  BICSQS  model  seem  reasonable.  The  second  and  third  modes  in  (28)  cause 
wave  speed  to  rise  slightly  at  higher  frequencies.  In  computing  the  modes  it 
was  assumed  that  in  equation  (28)  an  =  ci\/(2n  —  1)  and  kn  =  (2 n  —  l)2fci  with 
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Figure  18:  The  effect  on  compressional  wave  attenuation  of  varying  four  para¬ 
meters  in  the  BICSQS  model:  Ky  =  1  ±  0.2  GPa,  /*,  =  10  ±  5  kHz,  =  50  ±  25 
kHz,  v  =  0.15±0.1.  Green,  blue  and  red  curves  correspond  to  the  lower,  middle 
and  upper  values  respectively. 
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Figure  19:  Effect  upon  compressional  wave  speed  assuming  the  more  general 
equations  (28),  (29)  and  (30)  given  by  Biot  in  [8].  It  is  assumed  that  Ky  = 
10 9,fk  =  10 kHz,  =  50kHz.  In  the  upper  left  panel  s  =  0.75  and  green, 
blue  and  red  are  kq  =  10~9±1  in  ascending  order.  In  the  upper  right  panel 
kq  =  10-9  and  s  =  0.75  ±  0.25  with  the  same  color  ordering.  In  the  lower  right 
panel  kq  =  10-9,  s  =  0.75  and  green,  blue  and  red  correspond  to  using  1,2  and  3 
modes  in  equation  (28)  respectively.  In  the  lower  right  panel  k0  =  10~9,  s  =  1.0 
and  ba  =  105, 106, 10'  for  green,  blue  and  red  respectively. 

di  =  Wfc,  fci  =  Ky.  For  the  air  bubble  resonance  it  was  assumed  Co  =  (1  —  <j>)Kj 1 
and  k~l  =  (fiK^1  so  that  at  angular  frequency  uj  =  0  equation  (30)  gives  the 
Reuss  model  composite  compressibility 

c=(l  -fiKf  +  tK-1.  (41) 

The  gas  volume  fraction  was  taken  to  be  <j>  =  5  x  10~6.  This  sufficed  to  give 
agreement  with  the  measured  low  frequency  wave  speeds  shown  in  Figure  1.  The 
bulk  modulus  of  air  was  taken  to  be  Ka  =  0.142  MPa,  based  on  a  sound  speed 
of  343  m/s  and  a  density  of  1.21  kg/m3.  The  approximate  resonant  frequency 
was  taken  to  be  fa  =  25  kHz  and  the  mass  ma  in  (30)  was  computed  from 
ma  =  ka/(2nfa)2.  The  damping  parameter  ba  controls  how  rapidly  the  wave 
speed  changes  with  respect  to  frequency. 
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Figure  20:  Effect  upon  compressional  wave  attenuation  assuming  the  more  gen¬ 
eral  equations  (28),  (29)  and  (30)  given  by  Biot  in  [8].  It  is  assumed  that 
Ky  =  109,  fk  =  10  kHz,  /M  =  50 kHz.  In  the  upper  left  panel  s  =  0.75  and  green, 
blue  and  red  are  Kq  =  10~9±1  in  ascending  order.  In  the  upper  right  panel 
kq  =  10-9  and  s  =  0.75  ±  0.25  with  the  same  color  ordering.  In  the  lower  right 
panel  kq  =  10  9,  s  =  0.75  and  green,  blue  and  red  correspond  to  using  1,2  and  3 
modes  in  equation  (28)  respectively.  In  the  lower  right  panel  kq  =  10~9,  s  =  1.0 
and  ba  =  105, 106, 10'  for  green,  blue  and  red  respectively. 
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Sediment 

Pr 

Kr 

/? 

k 

a 

a 

Fine  Sand 

2670 

4.0  x  1010 

0.43 

3.12  x  10"14 

1.25 

1.20  x  10"6 

Medium  Sand 

2690 

3.2  x  101U 

0.385 

2.5  x  10'11 

1.35 

6.28  x  10-5 

Gravel 

2680 

4.0  x  1010 

0.30 

2.58  x  10"10 

1.25 

1.31  x  10"4 

Silty  Sand 

2670 

3.8  x  101U 

0.65 

6.33  x  10"ib 

3.0 

4.25  x  10"7 

Silty  Clay 

2680 

3.5  x  1010 

0.68 

5.2  x  10-14 

3.0 

1.24  x  10"6 

Table  2:  Biot-Stoll  parameters  for  five  different  sediment  types. 


Symbol 

Estimate 

PS 

1000 

Kf 

2.4  x  10y 

V 

1.01  x  10"3 

Table  3:  Typical  fluid  parameters  for  the  Biot-Stoll  model. 


2.6  Confirmation  of  Biot’s  equations 

The  Biot  model  is  often  described  as  phenomenological:  it  does  not  delve  into 
underlying  causes.  It  is  worth  noting  however  that  Burridge  and  Keller  [14] 
did  succeed  in  deriving  Biot’s  macroscopic  equations  from  the  microstructure, 
that  is,  using  the  linearized  equations  of  elasticity  and  fluid  motion  by  means  of 
the  mathematical  technique  of  two-scale  homogenization  under  the  assumptions 
that  the  medium  was  macroscopically  homogeneous  and  that  microstructure 
scale  was  much  smaller  than  the  macroscopic  scale.  They  also  demonstrated 
the  existence  of  a  strain  energy  function,  as  postulated  by  Biot. 


3  Determination  of  the  Biot-Stoll  parameters 

To  make  predictions  using  the  Biot-Stoll  model  the  thirteen  parameters  in  Ta¬ 
ble  1  must  be  determined.  Historically  this  has  been  done  eclectically  by  a 
combination  of  handbook  values,  empirical  formulas,  in  situ  measurements  and 
reference  to  the  literature.  The  exception  to  this  is  SAX99  in  which  attempts 
were  made  to  measure  most  of  the  Biot-Stoll  parameters  at  a  site  near  Fort 
Walton  Beach  Florida  (Table  1).  In  order  to  examine  the  predictions  of  the 
Biot  model  it  is  useful  to  have  parameters  for  a  variety  of  seabeds.  Values  for 
six  Biot-Stoll  parameters  for  five  different  sediments  are  given  in  Table  2.  The 
medium  sand  parameters  are  those  of  Table  1  for  the  sand  examined  in  SAX99. 
The  fine  sand,  silty  clay  and  gravel  parameters  are  from  Holland  and  Brunson 
[23].  The  parameters  for  silty  sand  are  from  Beebe,  McDaniels  and  Rubano  [2]. 
For  the  parameters  other  than  medium  sand  we  will  use  the  fluid  parameters  in 
Table  3. 

As  can  be  seen  in  Table  1  some  of  the  Biot-Stoll  parameters  can  be  deter¬ 
mined  quite  accurately  and  vary  little  in  the  context  of  water-saturated  seabeds. 
Let  us  examine  some  of  the  parameters  which  are  more  difficult  to  determine. 
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3.1  Porosity 

In  SAX99  porosity  was  estimated  using  several  different  techniques:  gravimetric 
determination  produced  a  95%  confidence  interval  [0.359,0.387].  Image  analysis 
produced  a  range  of  estimates  from  0.34  to  0.48.  Gravimetric  analysis  of  resin 
impregnated  cores  resulted  in  estimates  ranging  from  0.41  to  0.52.  Two  separate 
conductivity  probes  in  situ  yielded  95%  confidence  intervals  of  [0.364,0.383]  and 
[0.363,0.394].  These  and  all  subsequent  values  for  SAX99  are  from  Williams  et 
al.  [44].  Historically  porosity  has  been  measured  from  core  and  grab  samples. 
The  wide  range  of  values  generated  by  the  different  techniques  used  in  SAX99 
is  a  surprise  since  porosity  has  usually  not  been  regarded  as  being  subject  to 
great  uncertainty. 

3.2  Tortuosity 

For  SAX99  the  following  definition  of  tortuosity  was  used:  the  square  of  the  ratio 
of  the  minimum  length  of  a  contiguous  path  through  the  pore  space  between 
two  points  to  the  linear  distance  between  the  points.  Tortuosity  was  estimated 
by  image  analysis  of  resin-impregnated  divers  cores. 

Traditionally  experimenters  have  simply  followed  the  suggestion  of  Stoll  [37] 
and  used  a  =  1.25  for  highly  permeable  sediments  and  a  =  3.0  for  silty,  less 
permeable,  sediments  (see  Table  2). 

3.3  Fluid  bulk  modulus 

The  bulk  modulus  of  the  seawater  was  calculated  from 

Kf  =  pfc20 

where  c0  is  the  sound  speed  of  water  in  the  ocean.  The  density  p /  was  calculated 
based  on  the  temperature  and  salinity  of  the  water.  It  was  felt  that  assuming 
the  properties  of  the  pore  water  were  the  same  as  that  of  the  ocean  water  above 
introduced  little  error  ([44]).  The  result  was  the  rather  narrow  range  of  values 
shown  in  Table  1. 

Chotiros  and  Isakson  [17]  have  recently  challenged  this  range,  arguing  that 
the  effective  bulk  modulus  may  have  been  significantly  lower,  2.15  GPa  versus 

2.4  GPa,  due  to  the  presence  of  air  bubbles  in  the  pore  space.  As  evidence  for 
this  they  cite  the  attempts  to  measure  gas  content  of  sediment  cores  described 
in  Richardson  et  al.  [33]  which  resulted  in  measurements  ranging  from  20  to 
150  ppm.  The  authors  of  [33]  state  however  that  "All  of  these  values  indicate 
smaller  volumes  of  gas  than  the  system  was  designed  to  resolve.  We  cannot  on 
the  basis  of  these  data  reject  the  hypothesis  that  the  gas  content  is  zero."  As 
indicated  in  Figure  19  a  very  small  gas  fraction,  around  5  ppm,  would  suffice  to 
lower  the  low  frequency  wave  speed  the  measured  values  shown  in  Figure  1  if 
the  Reuss  model  (41)  is  assumed  for  the  effective  fluid  modulus5. 

5 Chotiros  and  Isakson  state  that  a  value  of  35.4  ppm  would  give  their  chosen  value  of 
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3.4  Permeability  and  the  pore  size  parameter 

For  SAX99  a  constant  head  technique  on  divers  cores  produced  a  95%  confidence 
interval  of  [2.1, 4.5]  x  10~n  m2  for  permeability.  Image  analysis  of  resin  impreg¬ 
nated  cores  produced  a  confidence  interval  of  [0.75,4.8]  x  1CT11  nr2.  An  in  situ 
constant  head  technique  produced  a  confidence  interval  of  [0.3,  6.1]  x  10-11  m2. 
Thus  the  estimates  spanned  more  than  a  decade.  As  mentioned  above  Williams 
et  al.  [44]  used  equation  (33)  to  determine  the  pore  size  parameter  from  the 
permeability  and  tortuosity. 

Traditionally  experimenters  have  calculated  permeability  from  statistics  on 
grain  size  distribution  via  a  variety  of  empirical  formulae.  For  instance  Holland 
and  Brunson  [23]  used  the  Kozeny-Carmen  equation 


k  = 


P3 

KS2{  l-/?)2 


(42) 


where  K  is  an  empirical  constant  which  is  approximately  5  for  spherical  grains 
and  So  is  the  surface  area  per  unit  volume  of  the  particles.  The  latter  parameter 
was  calculated  as 


so  =  E 

n 


6 

dn 


from  a  discrete  set  (dn,  wn)  of  grain  sizes  dn  and  proportions  wn  of  total  volume 
obtained  by  sorting  the  sample.  Beebe  et  al.  [2]  used  a  different  empirical 
relation  due  to  Krumbien  and  Monk  [27]  which  depends  upon  the  mean 
and  standard  deviation  cr^  of  grain  sizes  in  </>  units 


k  =  7.6 d2e-1-31<7*  x  10"lom2 


(43) 


where  d  =  2  M<t‘  mm  is  the  mean  grain  diameter.  Hovem  and  Ingram  [24] 
identified  the  pore  size  parameter  with  twice  the  hydraulic  radius  to  arrive  at 

“=3(nb)xlir3m-  <44) 

Permeability  is  then  calculated  from  the  Kozeny-Carmen  formula  written  as 


k  = 


pa2 
4 K' 


(45) 


Assuming  that  grain  sizes  have  a  ^-normal  distribution  (37)  Chotiros  [15]  arrives 
at  the  formula 

a  =  3(1 -/g)  exp(yo  +  x  10”3  m  (46) 

Kf  =  2.15  GPa,  but  do  not  state  what  model  this  is  based  upon.  This  may  be  a  misstatement 
since  the  rest  of  their  gas  fraction  values  are  roughly  consistent  with  the  Reuss-Voigt-Hill 
model  which  averages  the  Ruess  model  value  KJ1  =  (1  —  fyKw1  +  with  the  Voigt 

model  value  Kf  =  (1  —  4>)KW  +  (pKa  where  Kw  and  Ka  are  the  bulk  moduli  of  water  and  air 
respectively. 
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-12 
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0.85 

0.83 
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X 

io- 

-11 

2.4 

X 

io- 

-10 

4.7 

X 

io- 

-10 

1.3 

X 

io- 

-10 

0.38 

0.5 

1.0 

1.0 

X 

io- 

-io 

4.0 

X 

io- 

-io 

1.0 

X 

io- 

-y 

1.5 

X 

io- 

-io 

0.38 

1.2 

0.6 

6.6 

X 

io- 

-11 

1.5 

X 

io- 

-10 

2.1 

X 

io- 

-10 

1.1 

X 

io- 

-10 

Table  4:  Permeabilities  calculated  from  the  same  grain  size  statistics  using  four 
different  formulas. 


where  yo  =  —  in  2,  ay  =  —  in  2.  Finally  working  from  the  hydraulic  radius 
concept  Chotiros  [15]  arrived  at  a  formula  for  pore  size  which  differs  from  (46) 
only  by  a  minus  sign 

a  =  3(1  _  ~  exp(j/0  -  erg)  x  10~3  m .  (47) 

Table  4  shows  the  result  of  calculating  permeability  using  formulas  (43),  (44), 
(46),  and  (47)  (under  the  headings  Ukm,  km,  k^,  kc  respectively).  All  of  the 
grain  size  statistics  and  porosities  were  taken  from  [2],  Table  2.  As  can  be  seen 
the  four  formulas  lead  to  widely  differing  values  for  permeability,  especially  as 
the  grain  size  gets  smaller  and  the  dispersion  larger.  Thus  unfortunately  the 
value  of  this  very  important  parameter  is  elusive. 

3.5  Frame  response  parameters 

The  frame  response  parameters  are  the  shear  modulus  y*  and  the  frame  bulk 
modulus  K£.  They  are  usually  taken  to  be  complex  valued  in  order  to  account 
for  losses  due  to  intergranular  friction.  There  are  various  empirical  formulas  for 
calculating  them.  Bryan  and  Stoll  [10]  assume  a  functional  form  for  the  shear 
modulus 

yi  =  paa  exp  (-be)  (a0/pa)n  (48) 

where  pa  is  the  atmospheric  pressure, 


is  the  voids  ratio,  and  the  mean  effective  stress  due  to  over-burden  pressure  is 

=  1  +^K°  j  50“  P(z))(Pr~  Pf)dz  (49) 

where  z  is  the  depth  into  the  sediment,  K0  is  the  coefficient  of  earth  pressure  at 
rest,  which  is  typically  taken  to  be  0.5,  and  g  is  the  acceleration  due  to  gravity. 
Based  on  statistical  regressions  on  laboratory  results  Stoll  arrived  at  the  values 
a  =  2526,  b  =  1.504  and  n  =  0.448.  As  field  tests  tended  to  lead  to  somewhat 
higher  values  for  the  shear  modulus  than  laboratory  results,  Stoll  [38]  suggested 
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an  empirical  modification  /z  =  FF  ■  /zi  where  FF  =  2.  Badiey,  Cheng  and  Mu 
[1]  have  a  similar  formulation 

/z  =  a  x  105s~by/ao  (50) 

with  a  =  6.56,6  =  1.10  for  sand  dominant  sediments,  a  =  2.05,6  =  1.29  for 
clay  dominant  sediments  and  a  =  2.44,6  =  —1.628  for  general  (indeterminate) 
sediments.  The  complex  moduli  are  then  calculated  from 

li*  =  (1  +  iA^/n)  n  (51) 

where  Aj<b  and  are  the  log  decrements  for  compressional  and  shear  vibra¬ 
tions  and  v  is  the  Poisson  ratio. 

Figures  21-23  compare  the  predictions  for  the  variation  of  shear  wave  speed 
with  respect  to  depth  into  the  sediment  given  by  the  Bryan-Stoll  and  Badiey- 
Cheng-Mu  empirical  formulas  for  three  of  the  sediments  in  Table  2  for  which 
there  is  experimental  data.  The  data  points  for  Figure  21  are  from  [33]  and 
those  for  Figures  22  and  23  from  [34].  How  the  shear  wave  speeds  are  calculated 
will  be  discussed  in  Section  4.1.  It  should  be  noted  that  the  equations  of  motion 
were  derived  under  the  assumption  that  the  model  parameters  were  constant  and 
thus  use  the  equations  with  depth- varying  moduli  constitutes  an  approximation 
of  uncertain  accuracy.  Be  that  as  it  may,  Figures  21-23  indicate  that  the  Bryan- 
Stoll  formulation  with  either  FF  =  1  (silty  clay  and  fine  sand)  or  FF  =  2 
(medium  sand)  best  fits  the  data.  Possibly  the  Held  factor  is  needed  only  for 
highly  permeable  sediments. 

The  shear  wave  speed  is  significantly  affected  by  only  three  Biot-Stoll  para¬ 
meters:  the  shear  modulus  /z  =  Re/z*;  the  porosity  /?;  and  the  grain  density  pr 
[16].  In  SAX99  the  shear  wave  speed  was  measured  in  the  range  97  to  147 m/s 
and  the  average  value  cs  =  120  was  used  in  the  formula  p.  =  pc2s. 

4  Predictions  of  the  Biot-Stoll  model 

In  the  examples  of  the  predictions  of  the  Biot-Stoll  model  in  this  section  the 
frame  response  moduli  were  computed  from  (51)  at  a  depth  of  2  =  0.2  m  using 
the  empirical  formula  of  Bryan  and  Stoll  given  above.  The  required  Poisson 
ratios  and  log  decrements,  and  field  factors  are  given  in  Table  5. 

4.1  Wave  speed  and  attenuation 

Upon  taking  the  divergence  of  both  equations  in  (14)  and  assuming  time  har¬ 
monic  vibrations  e{x,  y,  z,  t)  =  e(x,  y,  z)e~lut, ...  we  obtain 

HV2e-CV2(  =  -uj2pe  +  uj2pf(  (52) 

CV2e  —  MV2C  =  — ui2pfe  +  (u>2m  +  iuid )  (. 
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Figure  21:  Predictions  for  the  variation  of  shear  wave  speed  with  respect  to 
depth  for  medium  sand  by  four  different  empirical  formulas  for  the  shear  wave 
modulus.  The  measured  values  are  from  [33].  The  curves  labeled  "BCM:  Gen¬ 
eral"  and  "BCM:  Sand"  were  calculated  from  (50)  using  the  parameters  given 
in  the  text.  The  curves  labeled  "BS:  FF  =  1"  and  "BS:  FF  =  2"  were  calculated 
from  (48)  with  the  indicated  field  factor. 


Sediment 

V 

FF 

Fine  sand 

0.25 

0.15 

0.15 

1 

Medium  sand 

0.15 

0.15 

0.3 

2 

Gravel 

0.25 

0.15 

0.15 

2 

Silty  sand 

0.25 

0.1 

0.1 

1 

Silty  clay 

0.25 

0.5 

0.5 

1 

Table  5:  Poisson  ratios  and  log  decrements  for  five  sediments. 
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Figure  22:  Predictions  for  the  variation  of  shear  wave  speed  with  respect  to 
depth  for  fine  sand  by  four  different  empirical  formulas  for  the  shear  wave  mod¬ 
ulus.  The  measured  values  are  from  [34],  The  curves  labeled  "BCM:  General" 
and  "BCM:  Sand"  were  calculated  from  (50)  using  the  parameters  given  in  the 
text.  The  curves  labeled  "BS:  FF  =  1"  and  "BS:  FF  =  2"  were  calculated  from 
(48)  with  the  indicated  field  factor. 


38 


Figure  23:  Predictions  for  the  variation  of  shear  wave  speed  with  respect  to 
depth  for  silty  clay  by  four  different  empirical  formulas  for  the  shear  wave  mod¬ 
ulus.  The  measured  values  are  from  [34],  The  curves  labeled  "BCM:  General" 
and  "BCM:  Clay"  were  calculated  from  (50)  using  the  parameters  given  in  the 
text.  The  curves  labeled  "BS:  FF  =  1"  and  "BS:  FF  =  2"  were  calculated  from 
(48)  with  the  indicated  field  factor. 
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Thus 

where 


I<  =  - 


H  -C 
C  -M 


— ui2  p  u2pf 

— cu2  Pf  oj2m  +  icod 


(53) 


Upon  taking  the  curl  of  both  equations  and  assuming  time-harmonic  vibrations 
we  obtain 


p\72ip  =  —  u2  (pip  —  pf&)  (54) 

0  =  —  ui2  (pf<p  —  m&)  +  dico& 

where  ip  =  V  X  u  and  =  V  X  w.  Solving  the  second  equation  for  and 
substituting  into  the  first  equation  gives 

V2p  +  k2p>  =  o 

where 

/  ui2(piom  +  ipd  —  p2ui ) 
s  y  (com  +  id)p 

Finally  solving  the  time-harmonic  version  of  (14)2  for  w  to  obtain 

w  =  — - - - —  [CVe  —  MV  C  +  u)2pf-u\ 

u>(u>m  +  id)  1  1  J 

and  substituting  this  result  into  (14)!  leads  to 

V2u  +  A1Ve  +  A2V(  +  k2u  =  0  (55) 

where 

^  ( H  —  pb)ojm  —  Cuipf  +  ( H  —  p)di 

1  ( urn  +  id)p 

A  MpfU)  —  Cnuu  —  iCd 
(com  +  id)p 

There  are  three  complex  wave  numbers  k\,  ^2,  ks  where  k\,  fc2  are  the  square 
roots  of  the  eigenvalues  of  K.  Thus  the  Biot-Stoll  model  predicts  two  com- 
pressional  waves  as  well  as  a  shear  wave.  The  two  compressional  waves  are 
referred  to  as  Type  I  and  II  or  fast  and  slow  waves.  Figures  24  and  25  show 
the  Biot-Stoll  models  predictions  for  Type  I  wave  speed  and  attenuation  for  the 
five  sediments  of  Table  2.  Type  I  waves  correspond  to  the  compressional  waves 
predicted  by  the  Buckingham  or  elastic  models,  but  for  highly  permeable  sedi¬ 
ments  such  as  the  gravel  and  SAX99  sediments,  they  are  more  dispersive  with 
respect  to  frequency.  The  behavior  of  the  attenuation  with  respect  to  frequency 
/  for  highly  permeable  sediments  is  roughly  characterized  as  increasing  like  / 2 
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Figure  24:  Biot  model  predictions  for  Type  I  compressional  wave  speeds  for  five 
sediments. 


for  low  frequencies  and  / 1/2  for  high  frequencies,  as  opposed  to  f1  in  the  Buck¬ 
ingham  model.  For  less  permeable  sediments  the  predictions  of  the  Biot  model 
are  similar  to  those  of  the  elastic  model.  Figures  26-29  show  the  Biot  model’s 
predictions  for  Type  II  and  shear  wave  speed  and  compressional  attenuation  for 
the  five  sediments  of  Table  2.  Both  waves  have  similar  speeds.  Type  II  waves 
for  unconsolidated  sediments  are  highly  attenuated  and  have  not  been  detected 
with  certainty.  See  however  Section  5.3. 

4.2  Reflection  and  transmission  at  the  ocean-sediment  in¬ 
terface 

4.2.1  The  case  of  constant  seabed  parameters 

The  equations  for  reflection  and  transmission  for  the  Biot  model  were  originally 
worked  out  by  Stoll  and  Kan  [40].  Stern,  Bedford  and  Milwater  [36]  have  con¬ 
sidered  the  case  of  depth- varying  sediment  properties.  We  will  treat  the  case  of 
constant  properties,  but  in  a  way  which  can  be  extended  in  theory  to  the  depth- 
varying  case.  Assuming  an  incident  plane  wave  of  amplitude  one,  pressure  in 
the  ocean  has  the  form 


Pa  =  eik*x  (eik*z  +  Re~ik‘z) 
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Figure  25:  Biot  model  predictions  for  Type  I  compressional  wave  attenuation 
for  five  sediments. 


where  kx  =  (w/co)  cos 9,  kz  =  (w/co)sin0,  cq  is  the  sound  speed  in  water,  the 
angle  of  incidence  6  is  measured  from  horizontal  and  the  positive  z-direction 
is  into  the  sediment.  The  time-harmonic  factor  e~lult  has  been  discarded  so 
that  P0  will  have  the  same  form  as  the  solutions  to  (53)  and  (55).  Vertical 
displacement  in  the  ocean  is  given  by 

i  ap  i 

Uzo  =  ——^  =  —  eikxX  (ikzeik‘z  -  Rikze~ik~z)  . 

pUJz  OZ  p(jJz 

For  the  ocean-sediment  interface  at  z  =  0,  the  following  conditions  are  imposed: 

Uzo(x,  0)  =  l3Uz(x,  0)  +  (1  -  /3)uz(x,  0)  =  uz(x,  0)  -  wz(x,  0)  (56) 

Po{x,  0)  =  tzz(x,  0) 

Po{x,  0)  =pf(x,  0) 

Txz{x,  0)  =  0. 

These  conditions  represent  continuity  of  displacement  in  the  ocean  with  aggre¬ 
gate  displacement  in  the  sediment,  continuity  of  normal  pressure,  continuity  of 
fluid  pressure,  and  vanishing  of  shear  stress  at  the  sediment  surface. 

Snell’s  law  requires  the  solutions  to  (53)  and  (55)  in  the  seabed  have  the 
form 


e(x,  z) 

_  0ikxx 

e(z) 

ux(x,z) 

ikxx 

ux{z) 

—  o 

—  fc- 

3 

42 


Figure  26:  Biot  model  predictions  for  Type  II  compressional  wave  speeds  for 
five  sediments. 


Substituting  the  first  of  these  into  (53)  gives 

0  ' 

0  ' 

For  a  constant  matrix  K  this  system  of  differential  equations  has  the  solution 

=  C\EieUlZ  +  C2E2eihz  (58) 

where  the  wave  numbers  are  given  by 

in  =  VK-  n  =  1,2 

and  \n,  En  are  eigenvalues  and  eigenvectors  of  K.  The  branch  cut  for  the  square 
root  is  chosen  so  that  Im  in  >  0.  The  exponentially  growing  components  of  the 
solution  (58)  involving  e~tlnZ  have  been  discarded.  Substituting  the  assumed 
form  for  u  from  (57)  into  (55)  gives 

u"  +  (k2s  -  k2x)u  =  -Ai  Ve  -  A2VC 

where  the  right  hand  side  is  known  by  virtue  of  (58).  This  system  has  solutions 
of  the  form 

Ux(z)  =  C3eu°z  +  FieUlZ  +  F2eU2Z,uz{z )  =  CAeu‘z  +  +Gie^lZ  +  G2eu*z 


e"{z) 

C(z) 


+  (K  —  k2I) 


e(z) 

C(z) 
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Figure  27:  Biot  model  predictions  for  Type  II  compressional  wave  attenuation 
for  five  sediments. 
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Figure  28:  Biot  model  predictions  for  shear  wave  speeds  for  five  sediments. 
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Figure  29:  Biot  model  predictions  for  shear  wave  attenuation  for  five  sediments. 


where  FA.F2.  G] ,  G2  are  linear  functions  of  Ci,C2  and  the  components  of  the 
eigenvectors  and  the  wave  number  is  given  by 

4  =  Vka  -  kl- 

By  equating  the  coefficients  of  elilZ,  e^2Z,  eltaZ,  the  constants  C3  and  C4 
can  be  related  and  the  components  of  the  eigenvectors  determined.  Thus  the 
four  sediments  solutions  e,(z),((z),ux(z),uz(z)  depend  upon  three  arbitrary 
constants.  The  four  interface  conditions  (56)  now  constitute  a  system  of  four 
linear  equations  which  may  be  solved  for  C\,  C2,  C3,  R. 

Figures  30  and  31  show  that  the  Biot  model  does  not  predict  a  critical 
angle  as  does  the  fluid  model,  but  it  does  predict  a  quasi-critical  angle  with  a 
similar  value.  For  instance  the  critical  angle  for  the  medium  sand  sediment  of 
SAX99  was  30°.  Figures  32  and  33  show  how  the  reflection  coefficient  varies 
with  frequency  at  an  incident  angle  of  5°.  Figure  34  shows  how  the  reflection 
coefficient  varies  with  frequency  at  different  angles  for  the  medium  sand  of  Table 
2. 

Given  that  the  Biot  model  predicts  that  transmission  of  energy  into  the 
seabed  is  possible  at  all  angles,  the  question  arises  as  to  whether  this  might 
explain  the  instances  of  detection  described  in  [41]  during  SAX99  of  objects 
by  waves  incident  at  angles  below  the  critical  angle.  Figure  31  indicates  that 
at  the  incident  angle  of  5°  and  frequency  20  kHz  used  in  the  experiment,  the 
coefficient  of  reflection  was  about  0.95.  Thus  other  explanations  such  as  diffrac- 
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Figure  30:  Biot  model’s  prediction  for  the  coefficient  of  reflection.  The  fre¬ 
quency  is  20  kHz  and  the  sound  speed  in  water  is  1460  nr  /  s,  which  is  below  the 
compressional  wave  speeds  of  all  of  the  sediments. 


tion  or  refraction  of  energy  into  the  sediment  by  surface  roughness  or  scattering 
of  evanescent  waves  by  volume  heterogeneity  within  the  sediment  seem  more 
plausible  [41],  however  see  Section  5.3. 


5  Difficulties  with  and  controversies  about  the 
Biot-Stoll  model 

5.1  Can  the  parameters  be  determined  accurately  enough? 

The  Biot  model  predictions  for  the  SAX99  data  shown  in  Figures  1  and  2  are 
best  fit  predictions  and  thus  only  can  be  obtained  a  posteriori.  It  is  worthwhile 
to  inquire  what  is  the  range  of  possible  predictions,  given  the  parameter  ranges 
in  Table  1,  Column  3.  Of  the  Biot-Stoll  parameters  for  which  there  is  much  un¬ 
certainty  only  four,  porosity,  permeability,  grain  bulk  modulus,  and  tortuosity, 
have  much  influence  on  Type  I  compressional  wave  speeds  and  attenuation6. 
Figures  35  and  36  show  the  influence  of  each  of  these  parameters  separately  on 

(> If  one  accepts  the  air  bubble  hypothesis  of  Chotiros  and  Isakson  [17],  then  fluid  bulk 
modulus  should  be  added  to  this  list.  See  Section  3.3  for  a  discussion  of  this  hypothesis  and 
Section  5.2.2  to  see  the  effect  of  varying  the  fluid  bulk  modulus. 
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Figure  31:  Biot  model’s  prediction  for  the  coefficient  of  reflection.  The  fre¬ 
quency  is  20  kHz  and  the  sound  speed  in  water  is  1530  m  /  s,  which  is  above  the 
compressional  wave  speeds  of  the  two  silty  sediments. 


47 


1.02 


1 


0.84 . . . . . . . . . . 

10'2  10'1  10°  10’  102  10" 
Frequency  (kHz) 


Figure  32:  The  reflection  coefficient  as  a  function  of  frequency  at  an  incident 
angle  of  5°.  The  sound  speed  in  water  was  1460m / s. 


compressional  wave  speed  and  attenuation  for  the  SAX99  sand. 

An  important  aspect  of  SAX99  was  that  most  of  the  parameters  were  mea¬ 
sured  sufficiently  many  times  to  permit  estimation  of  parameter  ranges.  Figures 
37  and  38  show  the  range  of  predictions  that  might  result  due  to  the  parameter 
ranges  given  in  Table  1.  In  particular  if  the  parameter  values  were  taken  to  be 
the  midpoints  of  the  ranges,  the  estimates  of  compressional  wave  speeds  would 
be  about  2%  too  high.  In  the  worst  case  in  which  the  data  lay  on  the  high  curve 
of  Figure  37,  the  midpoint  values  would  have  underestimated  the  wave  speed 
by  about  3.5%.  On  the  other  hand  the  predictions  for  wave  attenuation  would 
be  about  the  same  for  any  choice  of  values  within  the  range. 

5.2  The  frame  question 

There  are  concerns  about  the  Biot-Stoll  model,  insofar  as  its  applicability  to 
unconsolidated  seabeds,  centered  around  the  extent  to  which  an  uncemented 
collection  of  grains,  possibly  not  all  of  the  same  type  of  mineral,  can  be  charac¬ 
terized  as  an  elastic  frame. 

5.2.1  Is  the  model  applicable  to  heterogeneous  sediments? 

The  unjacketed  test  assumes  that  a  single  modulus  Kr  characterizes  the  frame 
material,  which  is  clearly  reasonable  only  for  a  homogeneous  frame,  that  is, 
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Figure  33:  The  reflection  coefficient  as  a  function  of  frequency  at  an  incident 
angle  of  5°.  The  sound  speed  in  water  was  1530m / s. 


one  composed  of  grains  of  a  single  mineral.  Biot  and  Willis  [5]  conjectured 
that  a  heterogeneous  material  might  act  as  an  "equivalent  homogeneous  solid", 
however  Hickey  and  Sabatier  [22]  point  out  this  would  require  that  changes  in 
grain  shape  due  to  deviatoric  strains  cancel  out,  which  may  be  optimistic.  Thus 
the  question  of  what  to  do  about  heterogeneous  materials  remains  open. 

5.2.2  Is  porosity  constant? 

Another  question  is  whether  in  the  unjacketed  test  porosity  is  unchanged  when 
the  specimen  is  subjected  to  fluid  pressure.  Quite  possibly  for  an  uncemented 
collection  of  grains  it  would  change.  In  this  case  the  moduli  in  (8)  would  be 
calculated  using  (9)  rather  than  (11),  but  whereas  the  bulk  modulus  Kf  of  water 
is  well  known,  it  is  not  clear  how  to  measure  the  modulus  Kfr.  Inversions  on 
four  data  sets  by  Chotiros  [16]  gave  values  of  Kfr  that  were  60%,  62%,  88%,  61% 
of  the  values  Kjr that  would  be  obtained  by  use  of  (10).  The  88%  figure  was 
from  data  on  uncompacted  laboratory  sand.  He  calls  the  hypothesis  Kfr  ^  K®r 
the  independent  coefficient  of  fluid  content  hypothesis  (ICFC).  As  is  indicated 
in  Figure  1,  the  Biot  model  prediction  for  compressional  wave  speed  was  higher 
than  the  measured  values  for  the  two  measurements  below  one  kHz.  To  see  if 
the  ICFC  hypothesis  might  explain  this  let  Kfr  =  6K ®r.  Figures  39  and  40 
show  the  results  of  varying  9  on  compressional  wave  speed  and  attenuation. 


49 


Frequency  (kHz) 

Figure  34:  Reflection  coefficient  as  a  function  of  frequency  for  different  incident 
angles  for  medium  sand.  The  critical  angle  is  around  30  ° . 


A  value  of  around  9  =  0.95  sufficed  to  lower  the  low  frequency  predictions  to 
near  the  measured  values,  but  reduced  the  predictions  at  higher  frequencies 
to  well  below  the  measured  values.  However  as  Figures  41  and  42  show  this 
can  be  ameliorated  by  increasing  the  permeability  from  k  =  2.5  x  10~n  m2  to 
k  =  5  x  10~n  nr2,  and  the  tortuosity  from  a  =  1.35  to  a  =  1.12  as  was  done 
to  generate  the  "best  fit"  plot  in  Figure  1.  Thus  this  hypothesis  seems  to  have 
merit.  It  provides  an  alternative  to  simply  raising  the  porosity  to  get  a  better 
fit  to  the  low  frequency  data.  On  the  other  hand,  while  there  is  a  hypothetical 
experiment  by  which  the  additional  parameter  Kfr  might  be  measured  in  Biot 
and  Willis  [5],  it  is  not  clear  whether  it  is  practical.  It  should  also  be  noted 
that  the  air  bubble  hypothesis  of  Chotiros  and  Isakson  [17]  (Section  3.3)  leads 
to  the  same  conclusion.  Indeed  if  we  accept  (10),  but  use  (41)  to  compute  Kf, 
then  the  resulting  value  of  D  computed  from  (11)  can  be  equated  to  that  of  (9), 
resulting  in  a  one-to-one  correspondence  between  the  gas  volume  fraction  and 
9.  Gas  fractions  of  0=1,3  and  5  ppm  correspond  to  values  of  9  =  0.984, 0.950 
and  0.919. 

5.2.3  Should  all  of  the  grains  be  included  in  the  frame? 

A  final  difficulty  with  treating  an  unconsolidated  sediment  as  an  elastic  frame  is 
the  problem  of  what  should  be  included  in  the  frame.  In  [28]  it  is  shown  that  the 
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Porosity  Permeability 


Grain  bulk  modulus  Tortuosity 


Figure  35:  Influence  of  four  Biot-Stoll  parameters  on  compressional  wave  speed. 
The  green  curve  used  the  lower  bound  for  the  parameter  given  in  Table  1, 
Column  3,  the  red  curve  the  upper  bound,  and  the  blue  curve  the  midpoint. 
The  midpoint  values  were  used  for  all  parameters  other  than  the  one  being 
varied. 
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Figure  36:  Influence  of  four  Biot-Stoll  parameters  on  compressional  wave  at¬ 
tenuation.  The  green  curve  used  the  lower  bound  for  the  parameter  given  in 
Table  1,  Column  3,  the  red  curve  the  upper  bound,  and  the  blue  curve  the 
midpoint.  The  midpoint  values  were  used  for  all  parameters  other  than  the  one 
being  varied. 
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Figure  37:  Range  of  predictions  for  compressional  wave  speed  possible  given  the 
parameter  ranges  of  Table  1,  Column  3.  The  predictions  labelled  "Best  in  box" 
are  those  for  the  parameters  of  Column  4.  The  predictions  labelled  "Best  out  of 
box"  are  those  of  Column  4  with  the  changes  (3  =  41.5,  a  =  1.12,  k  =  5  x  10~n. 
The  "Midpoint"  curve  used  the  midpoints  of  the  ranges  in  Column  3.  The 
"High"  and  "Low"  curves  were  obtained  by  moving  each  the  parameters  porisity, 
permeablity,  grain  bulk  modulus,  and  tortuosity  to  whichever  range  endpoint 
increased  or  decreased  the  wave  speed  respectively. 
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Figure  38:  Range  of  predictions  for  compressional  wave  attenuation  possible 
given  the  parameter  ranges  of  Table  1,  Column  3.  The  predictions  labelled  "Best 
in  box"  are  those  for  the  parameters  of  Column  4.  The  predictions  labelled  "Best 
out  of  box"  are  those  of  Column  4  with  the  changes  /3  =  41.5,  a  =  1.12,  k  = 
5  x  ICC11.  The  "Midpoint"  curve  used  the  midpoints  of  the  ranges  in  Column 
3.  The  "High"  and  "Low"  curves  were  obtained  by  moving  each  the  parameters 
porisity,  permeablity,  grain  bulk  modulus,  and  tortuosity  to  whichever  range 
endpoint  increased  or  decreased  the  wave  speed  respectively. 
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Figure  39:  Result  on  compressional  wave  speed  of  varying  the  parameter  6  in 
the  equation  Kfr  =  8K^r. 


55 


Figure  40:  Result  on  compressional  attenuation  of  varying  the  parameter  8  in 
the  equation  Kfr  =  6K^r. 
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Figure  41:  Result  on  compressional  wave  speed  of  varying  permeability  when 
Kfr  =  0K®r  with  9  =  0.95.  "bob"  is  the  best  fit  from  Figure  1. 
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Figure  42:  Result  on  compressional  wave  attenuation  of  varying  permeability 
when  Kfr  =  0Kjr  with  9  =  0.95.  "bob"  is  the  best  fit  from  Figure  2. 
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application  of  an  external  force  to  a  granular  lattice  causes  the  formation  of  force 
chains  (see  Figure  43).  Thus  it  is  not  clear  that  all  of  the  sedimentary  grains 
should  be  considered  part  of  the  frame.  Conversely  it  may  be  that  some  fluid  is 
trapped  in  cracks  in  the  grains  or  between  grains  and  should  be  considered  part 
of  the  frame.  Chotiros  [16]  defined  three  additional  porosities  pf,ps  and  pc. 
The  first  two  are  the  portion  of  fluid  trapped  in  the  frame  and  the  portion  of 
the  grains  to  be  regarded  as  part  of  the  fluid.  The  third,  the  composite  porosity, 
is  defined  by 


Pc 


P-Pf 

1  -Ps-Pf' 


(59) 


The  composite  densities  and  bulk  moduli  are  defined  to  be 


Prc 

K 

xvrc 


(1  -  fif)pr  +  fifPf,  pfc  =  (1  -  Ps)pf  +  Ppr 

KrKf  =  KrKf 

(1  -pf)Kf  +  pfKr'  sc  (1  —  Ps)Kr  +  PsKf ' 


(60) 


Inversions  on  four  data  sets  gave  values  of  pf  =  0.15, 0,  0,  0  and  /3S  =  0.36,  0.36, 
0.07,0.31  for  seabeds  of  porosity  p  =  0.37,0.36,0.41,0.40  giving  composite 
porosities  Pc  =  0.44,0.57,0.44,0.58  which  are  considerably  higher  than  the 
measured  porosities.  This  hypothesis,  which  Chotiros  calls  the  composite  ma¬ 
terials  (CM)  hypothesis,  did  not  seem  to  explain  the  higher  apparent  porosity, 
P  =  0.415,  for  the  SAX99  data.  However  since  the  bulk  moduli  Kf  and  Kr  oc¬ 
cur  only  in  the  unjacketed  test  (see  Section  2.1)  where  the  presence  or  absence 
of  force  chains  would  be  irrelevant,  it  is  not  clear  that  the  use  of  the  composite 
moduli  (60)2  is  required.  If  one  leaves  Kf  and  Kr  as  they  are  and  replaces 
all  instances  of  p,pj  and  pr  by  their  composites  given  in  (59)  and  (60)  1,  then 
the  resulting  predictions  are  shown  in  Figures  44  and  45.  It  should  be  noted 
that  these  graphs  are  sufficiently  similar  to  those  of  Figures  39  and  40  that  the 
composite  materials  and  independent  coefficient  of  fluid  content  hypotheses  are 
indistinguishable  on  the  basis  of  their  predictions  of  wave  speed  and  attenuation 
alone. 


5.3  The  "fast  slow  wave"  controversy 

When  in  three  experiments  waves  having  an  apparent  speed  of  about  1200  m  /  s 
were  detected  at  subcritical  grazing  angles,  Chotiros  [15]  conjectured  that  they 
were  Type  II  compressional  waves.  As  indicated  in  Figure  26  Type  II  waves 
are  not  predicted  to  have  speeds  nearly  that  large  with  the  parameter  values 
conventionally  used  in  the  Biot-Stoll  model.  To  arrive  at  a  Type  II  wave  speed 
of  1200  m/s  Chotiros  suggested,  based  upon  measurements  of  the  grain  bulk 
modulus  by  Molis  and  Chotiros  [31],  that  the  value  of  the  grain  bulk  modulus 
Kr  should  be  around  7  x  109  Pa  as  opposed  to  the  conventional  value  of  about 
4  x  1010  Pa,  which  is  the  bulk  modulus  of  quartz  crystals'  .  As  indicated  in 

'Chotiros  suggests  that  the  quartz  crystal  value  is  not  applicable  because  of  cracks  and 
other  imperfections  in  the  sand  grains. 
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Fig.  1  .  An  image  of  tire  force  chains  in  a  granular 
medium  as  viewed  between  two  crossed  circular 
polarlzors.  The  particles.  3  mm  Pyrex  spheres, 
are  surrounded  by  an  Index -matching  fluid,  a 
mixture  of  glycerol  and  water.  The  beads  are 
placod  in  a  box  with  aidos  70  mrn  by  70  mm  by 
40  mm;  the  image  size  is  40  mm  by  42  mrn.  A 
2O0-N  force  is  exerted  on  a  piston  that  covers 
the  top  surface .  The  stressed  boods  can  bo  soon 
as  bright  regions  (the  brightness  is  due  to  their 
stress -induced  birefringence). 
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Figure  43:  Illustration  of  a  force  chain  in  sand. 
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Figure  44:  Compressional  wave  speeds  for  the  SAX99  sand  when  the  proportion 
Ps  of  grains  to  be  included  in  the  pore  fluid  is  varied.  The  proportion  /3f  of  fluid 
included  in  the  frame  was  assumed  to  be  zero.  The  parameters  used  were  those 
of  Table  1,  Column  4. 
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Figure  45:  Compressional  wave  attenuation  for  the  SAX99  sand  when  the  pro¬ 
portion  j3s  of  grains  to  be  included  in  the  pore  fluid  is  varied.  The  proportion 
/3 f  of  fluid  included  in  the  frame  was  assumed  to  be  zero.  The  parameters  used 
were  those  of  Table  1,  Column  4. 
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Figure  46:  The  effect  of  using  the  non-conventional  parameter  values  Kr  = 
1  x  109  Pa,  Re  Kb  =  5.3  x  109Pa  suggested  by  Chotiros  on  the  predictions  for 
Type  I  waves  of  the  Biot  model  for  the  SAX99  sediment.  Blue  lines  are  the 
predictions  for  the  SAX99  best  fit  parameters  of  Table  1.  Red  lines  are  the 
predictions  resulting  from  the  alternative  parameters. 


Figure  35,  lowering  the  value  of  Kr  lowers  the  predicted  Type  I  wave  speeds. 
To  achieve  agreement  with  the  measured  Type  I  wave  speeds  in  the  three  ex¬ 
periments  Chotiros  assumed  a  frame  bulk  modulus  of  Kf,  =  5.3  x  109  Pa,  which 
is  about  50  times  larger  than  the  conventional  values.  Figures  46  and  47  show 
the  effects  of  these  modifications  on  the  predictions  for  the  SAX99  sand.  Type 
I  compressional  waves  become  less  dispersive  with  respect  to  frequency  and  the 
growth  of  attenuation  more  nearly  linear  with  respect  to  frequency.  Type  II 
waves  now  have  speeds  in  the  1200-1400  m/s  range  at  frequencies  above  1  kHz 
and  are  less  attenuated. 

Subsequent  work,  both  experimental  and  theoretical,  has  tended  to  refute 
Chotiros’  conjecture.  As  indicated  above  the  SAX99  measurements  confirmed 
the  conventional  values  (Table  1).  Simpson  and  Houston  [35]  found  in  laboratory 
experiments  that  for  a  smoothed  water-sediment  interface  no  waves  other  than 
those  having  the  expected  speeds  of  Type  I  waves  could  be  detected,  but  for  a 
roughened  interface  and  shallow  grazing  angles  a  wave  with  an  apparent  speed  of 
1200  m/s  was  detected.  The  authors  state  "However,  for  near  normal  incidence, 
these  measurements  have  enough  temporal  resolution  to  clearly  show  that  this 
later  time  arrival  is  better  described  by  the  roughened  interface  scattering,  which 
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Figure  47:  The  effect  of  using  the  non-conventional  parameter  values  Kr  = 
7  x  109Pa,  ReKb  =  5.3  x  109Pa  suggsted  by  Chotiros  on  the  predictions  for 
Type  II  waves  of  the  Biot  model  for  the  SAX99  sediment.  Blue  lines  are  the 
predictions  for  the  SAX99  best  fit  parameters  of  Table  1.  Red  lines  are  the 
predictions  resulting  from  the  alternative  parameters. 
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give  the  appearance  of  a  wave  front  with  a  virtual  wave  speed  of  1200  m/s  at 
shallow  ensonification  angles,  rather  than  by  a  1200  m/s  slow  compressional 
wave."  On  the  theoretical  side  Hickey  and  Sabatier  [22]  argue  on  the  basis 
of  the  Voigt-Reuss-Hill  formulae  that  for  the  mineral  composition  of  the  sand 
used  in  one  of  Chotiros’  experiments  the  grain  bulk  modulus  must  lie  between 

3.3  x  1010  and  3.7  x  1010  Pa.  This  range  does  not  encompass  Chotiros’  value  of 
Kr  =  7  x  109  Pa.  They  also  argue  that  the  constraint 

0  <1  —  Kb/Kr  <  1 

([5])  requires  that  Kb  <  4.2  x  109  Pa  which  again  is  not  consistent  with  Chotiros’ 
value  of  Kb  =  5.3  x  109  Pa. 

5.4  The  Biot  model’s  predictions 

5.4.1  The  controversy  about  the  growth  of  attenuation  with  increas¬ 
ing  frequency 

It  is  generally  agreed  that  internal  dissipation  generates  wave  attenuation  that 
is  linear  in  frequency  [26] .  The  Biot  model’s  prediction  that  compressional  wave 
attenuation  is  not  necessarily  a  linear  function  of  frequency,  but  is  more  like  f2 
at  low  frequencies  and  f1!2  at  high  frequencies  is  due  to  viscous  losses.  This 
prediction  has  been  a  source  of  controversy  for  nearly  three  decades.  Hamilton 
[21]  stated  "The  number  studies  of  attenuation  in  sands  from  1-100  kHz  affirm 
an  approximate  first-power  relationship  between  attenuation  and  frequency,  and 
deny  any  f1/2  relationship  as  called  for  by  Stoll  in  this  frequency  range."  More 
recently  Buckingham  [11]  asserted  "A  substantial  body  of  data  supports  the 
view  that  compressional-wave  attenuation  in  many  porous,  granular  materi¬ 
als  varies  more  or  less  accurately  as  the  first  power  of  frequency,  /  ,  over  an 
extended  frequency  range,  from  1  Hz  up  to  1  MHz."  On  the  other  hand  Kib- 
blewhite  stated  in  a  survey  article  [26]  "...  in  spite  of  the  sparsity  of  data  for 
marine  sediments  at  low  frequencies,  the  widely  held  assumption  that  the  rela¬ 
tionship  between  attenuation  and  frequency  is  linear  from  seismic  to  ultrasonic 
frequencies  does  not  appear  justified."  Brunson  and  Johnson  [9]  measured  shear 
wave  attenuation  in  the  1-10  kHz  range  and  found  "...  for  a  sediment  of  this 
type,  the  relative  motion  of  the  pore  fluid,  and  the  sediment  skeletal  frame  is  a 
significant  contributor  to  the  observed  loss  of  energy...".  Turgut  and  Yamamoto 
[43]  found  that  "Consistency  between  predicted  and  directly  measured  porosity 
and  permeability  indicates  that  the  viscous  loss  mechanism  has  a  major  effect 
on  attenuation  in  the  frequency  range  of  1-30  kHz  for  sandy  sediments."  The 
source  of  the  difficulty  is  illustrated  in  Figure  2.  It  not  easy  to  differentiate 
between  the  predictions  of  the  Biot  and  Buckingham  model  in  the  frequency 
range  1-100  kHz,  and  this  historically  has  been  the  range  in  which  most  mea¬ 
surements  of  attenuation  have  been  made.  The  SAX99  measurements  do  seem 
to  indicate  that  the  Biot  model’s  prediction  of  attenuation  varying  as  f1^2  at 
high  frequencies  is  wrong,  and  this  is  progress. 
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As  pointed  out  in  [44]  experimental  evidence  for  strong  velocity  dispersion,  as 
predicted  by  the  Biot  model,  is  also  evidence  that  attenuation  cannot  be  a  linear 
function  of  frequency.  Causality  dictates  the  Kramers-Kronig  integral  relations 
between  wave  speed  and  attenuation.  In  O’Donnel  et  al.  [32]  it  is  argued  that 
the  Kramers-Kronig  integral  relations  lead  to  the  following  approximate  local 
relation  between  wave  speed  c  and  attenuation  a  (in  nepers/m) 
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Thus  for  the  Buckingham  model  which  predicts,  to  the  first  order  in  w,  that 
attenuation  is  given  by  a  =  qc o,  wave  speeds  must  be  logarithmically  dispersive 
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Figure  48  shows  the  Kramers-Kronig  compatible  fit  to  the  SAX99  data  from 
[44] 
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with  Ci  =  1705  nr  /s,  uq/27r  =  100  Hz,  7  =  5.50  x  10_6Np-s  /m ,  r  =  0.3  ms,  So  = 
0.05  ms  /  nr,  as  well  as  the  prediction  of  the  Buckingham  model.  As  can  be  seen 
Kranrers-Kronig  relation  implies  that  a  strong  velocity  dispersion  necessitates  a 
departure  from  attenuation  that  is  linear  in  frequency.  There  is  other  evidence 
for  strong  velocity  dispersion.  Maguer  et  al.  [30]  have  found  that  experiments  on 
subcritical  penetration  of  energy  into  the  seabed  suggest  a  substantially  lower 
conrpressional  speed  than  is  measured  from  cores  at  200  kHz  and  that  these 
lower  speeds  are  compatible  with  those  predicted  by  the  Biot  model.  Turgut 
and  Yamamoto  [43]  also  found  velocity  dispersion  compatible  with  Biot  model 
predictions  in  cross-hole  tomography  experiments  (see  Figure  53).  Thus  while 
the  SAX99  data  tends  to  affirm  / 1  growth  at  high  frequencies,  it  adds  to  the 
evidence  that  low  frequency  attenuation  is  not  proportional  to  /  ' . 


6  Incorporation  of  frequency-dependent  viscoelas¬ 
tic  mechanisms 

As  indicated  in  Figures  1  and  2  the  "standard"  Biot  model  with  frequency- 
independent  bulk  and  shear  moduli  fails  to  agree  with  the  SAX99  data  in  two 
significant  ways.  It  predicts  wave  speeds  that  were  higher  than  measured  at  low 
frequencies,  at  least  if  the  parameters  are  confined  to  the  ranges  given  in  Table 
1,  Column  3.,  and  it  predicts  attenuations  that  are  too  low  at  high  frequencies. 
Ways  of  addressing  the  low  frequency  problem  have  been  addressed  above.  One 
can  posit  air  bubbles  in  the  pore  space  as  Chotiros  and  Isakson  did  in  [17] 
(Sections  2. 4, 3. 3).  Alternatively  one  can  adopt  the  independent  coefficient  of 
fluid  content  (ICFC)  hypothesis  discussed  in  Section  5.2.2  or  the  composite 
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1.2 


Figure  48:  Blue  line:  Kramers-Kronig  compatible  fit  to  SAX99  compressional 
wave  and  attenuation  data.  Red  line:  Predictions  of  the  Buckingham  model. 


material  (CM)  hypothesis  of  Section  5.2.3.  As  can  be  seen  in  Figures  44  and  39 
any  of  these  conjectured  mechanisms  have  the  effect  of  simply  translating  the 
wave  speed  curves  downward  and  thus,  absent  any  means  of  determining  the 
controlling  parameters  other  than  choosing  them  to  agree  with  the  measured 
data,  there  is  no  way  to  distinguish  among  them.  The  problem  of  the  poor 
predictions  for  attenuation  at  high  frequencies  has  been  addressed  recently  by 
means  of  the  BICSQS  model  discussed  in  Section  2.4.  In  this  section  we  examine 
the  extent  to  which  it  rectifies  this  problem. 

As  a  starting  point  we  took  the  midpoints  of  the  intervals  in  Table  1,  Col¬ 
umn  3  as  the  values  of  the  Biot-Stoll  parameter  (see  Figure  37).  The  Poisson 
ratio,  which  has  negligible  effect,  was  set  to  v  =  0.15,  following  [17].  The  shear 
modulus  was  adjusted  to  give  the  experimentally  measured  value  of  120  m/s 
[44],  The  Reuss  model  of  composite  compressibility  (41)  was  used  to  obtain 
agreement  with  the  low  frequency  wave  speed  data  by  specifying  a  gas  vol¬ 
ume  fraction  q i.  Agreement  with  the  high  frequency  measurements  for  wave 
speed  and  attenuation  were  obtained  by  manipulating  the  parameters  Ky  and 
ffi  respectively.  Finally  agreement  with  the  measured  values  of  wave  speed  and 
attenuation  in  the  1  —  10  kHz  range  was  obtained,  to  the  extent  possible,  by 
adjusting  the  bulk  relaxation  frequency  // .  A  compromise  was  required  here 
because  values  that  improved  agreement  with  attenuation  worsened  it  with  re¬ 
spect  to  wave  speed.  Figures  49  and  50  show  the  results  of  these  machinations 
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Parameter 

Units 

SAX-CI 

S  AX-Mid/Dps 

TY 

Porosity  f3 

0.37 

0.3785 

0.44 

Grain  density  pr 

kg/nr3 

2690 

2690 

2650 

Fluid  density  pf 

kg/rn3 

1023 

1023 

1000 

Grain  bulk  modulus  Kr 

GPa 

36 

40.5 

36 

Effective  fluid  bulk  modulus  Kf 

GPa 

2.15 

2.19 

2.25 

Fluid  viscosity  p 

kg/rn-s 

0.001 

0.00105 

0.001 

Permeability  k 

115 

33 

17.5 

Pore  size  a 

57 

31 

28 

Structure  factor  a 

1.35 

1.38 

1.24 

Frame  shear  modulus  p 

GPa 

0.028 

0.03 

0.024 

Gap  modulus  Ky 

GPa 

0.968 

0.8/0. 9 

0.9 

Bulk  relaxation  frequency  fk 

kHz 

4.8 

4/1.5 

3.2 

Shear  relaxation  frequency 

kHz 

56 

70/80 

“TO5 

Gas  fraction  <f> 

xlO~° 

6.7* 

5.5 

1.4* 

Pore  size  std  dev  a 

0 

0/1.25 

0 

Table  6:  SAX-CI:  Parameters  for  the  SAX99  sediment  used  by  Chotiros  and 
Isakson.  SAX-Mid:  Midpoints  of  the  parameter  ranges  in  Column  3  of  Table 
1.  The  alternative  Dps  values  in  Column  4  assume  distributed  pore  sizes.  TY: 
The  parameters  Chotiros  and  Isakson  used  to  model  the  cross-hole  tomogra¬ 
phy  experiemnt  of  Tugut  and  Yamamto.  The  parameters  marked  with  *  were 
calculated  using  the  Reuss  model  for  composite  compressibilities.  Chotiros  and 
Isakson  used  a  different  model.  The  TY  data  was  not  measured  at  sufficiently 
high  frequencies  to  ascertain  the  shear  relaxation  frequency.  Chotiros  and  Isak¬ 
son  listed  it  as  "infinity". 


along  with  the  curves  that  Chotiros  and  Isakson  arrived  at  by  using  the  same 
procedure.  The  fit  that  Chotiros  and  Isakson  arrived  at  is  somewhat  better, 
but  as  indicated  in  Table  6,  they  used  a  value  of  permeability  that  was  more 
than  three  times  the  midpoint  of  the  estimated  range  given  in  Table  1.  The 
only  rationale  they  give  for  this  adjustment  is  that  "This  is  consistent  with  ob¬ 
servations  that  permeability  measurements  may  only  be  accurate  to  within  an 
order  of  magnitude."  As  pointed  out  in  Section  3.4  permeability  was  measured 
in  three  different  ways  during  the  SAX99  exercise  and  none  of  them  arrived  at  a 
value  nearly  as  high  as  the  Clrotiros-Isakson  value.  This  being  the  case,  a  better 
alternative  may  be  to  use  the  distributed  pore  size  viscosity  correction  factor  of 
Yamamoto  and  Turgut  [45]  presented  in  Section  2.5.2.  Figures  51  and  52  show 
that  when  a  standard  deviation  of  a  =  1.25  is  assumed  and  the  mean  pore  size 
distribution  is  computed  from  (39)  and  the  midpoint  of  the  permeability  range 
in  Table  1,  agreement  with  the  SAX99  data  is  comparable  to  that  obtained  by 
Chotiros  and  Isakson. 

In  [17]  Chotiros  and  Isakson  discuss  the  predictions  of  their  BICSQS  model 
for  five  sets  of  data  taken  from  the  literature.  Other  than  the  SAX99  data 
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1.2 


1.18 
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Figure  49:  Best  fits  to  SAX99  wave  speed  data. 


Figure  50:  Best  fit  to  the  SAX99  attenuation  data. 


69 


1.2 


1.18 


Figure  51:  Best  fit  to  the  SAX99  wave  speed  data,  assuming  distributed  pore 
sizes. 


the  only  one  of  the  five  that  has  measurements  of  compressional  wave  speed 
against  frequency  are  the  cross-hole  tomography  experiments  of  Turgut  and 
Yamamoto  [43].  Figures  53  and  54  show  Chotiros  and  Isakson’s  best  fit  to  the 
wave  and  attenuation  data.  Table  6,  Col.  5  gives  the  parameters  they  used. 

It  may  be  noted  that  the  experimental  values  of  wave  speed  in  Figure  53 
ascend  with  frequency  more  quickly  than  is  predicted  by  the  BICSQS  model. 
Figure  55  shows  that  a  better  fit  to  the  wave  speed  data  can  be  obtained  by 
incorporating  Biot’s  model  for  air  bubble  resonance  (equation  (30))  discussed  in 
Section  2.4.  Figure  56  shows  that  with  or  without  bubble  resonance  the  BICSQS 
model  overestimated  the  compressional  wave  attenuation  in  the  lower  frequency 
range.  However  Turgut  and  Yamamoto  ([43],  Figure  6)  indicate  that  calculating 
the  wave  speeds  from  the  measured  attenuation  using  the  approximate  local 
Kramers-Kronig  relations  (61)  somewhat  underestimated  the  measured  wave 
speeds.  This  is  illustrated  in  Figure  57  where  it  is  seen  that  the  wave  speeds 
computed  from  (61)  using  numerical  integration  on  the  attenuation  data  in 
Figure  56  are  not  in  as  good  agreement  with  the  measured  wave  speed  values 
as  the  wave  speeds  calculated  from  (61)  using  the  higher  attenuation  values 
predicted  by  the  four-mode  bubble  resonance  model. 


70 


Figure  52:  Best  fit  to  the  SAX99  attenuation  data,  assuming  distributed  pore 
sizes. 


Parameter 

Units 

Cl 

Bubble  resonance 

Gap  modulus  Ky 

GPa 

1.35 

0.55 

Bulk  relaxation  frequency  fk 

kHz 

4.75 

3.2 

Shear  relaxation  frequency 

kHz 

10B 

100 

Gas  fraction  <j> 

xlO-6 

3.75 

3.75 

Pore  size  std  dev  a 

0 

0 

Resonant  frequency  fa 

kHz 

4.5 

Damping  parameter  ba 

MPa-s 

0.7 

Table  7:  Parameters  used  to  compare  the  BICSQS  model  with  (Col.  4)  and 
without  (Col  3)  air  bubble  resonance.  The  TY  data  was  not  measured  at  suf¬ 
ficiently  high  frequencies  to  ascertain  the  shear  relaxation  frequency.  Chotiros 
and  Isakson  listed  it  as  "infinity". 


71 


Figure  53:  Best  fit,  as  determined  by  Chotiros  and  Isakson  [17],  to  wave  speed 
data  of  the  cross- hole  tomography  of  Turgut  and  Yamamoto  [43].  Also  shown 
are  the  predictions  of  the  standard  Biot  model  with  =  43.6+2.08*  MPa,  /. i *  = 
29.2  +  1.8*  MPa. 
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Figure  54:  Best  fit,  as  determined  by  Chotiros  and  Isakson  [17],  to  wave  at¬ 
tenuation  data  of  the  cross-hole  tomography  of  Turgut  and  Yamamoto  [43]. 
Also  shown  are  the  predictions  of  the  standard  Biot  model  with  =  43.6  + 
2.08*  MPa,  n*  =  29.2  +  1.8*  MPa . 
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Figure  55:  The  predictions  of  the  BICSQS  model  when  air  bubble  resonance 
is  included.  The  parameters  used  are  given  in  Table  7.  The  curve  labeled 
Chotiros-Isakson  was  obtained  by  reducing  the  effective  fluid  bulk  modulus  to 
agree  with  the  value  which  gave  good  agreement  at  low  frequencies  in  the  bubble 
resonant  case  and  the  adjusting  the  gap  modulus  and  bulk  relaxation  frequencies 
to  improve  the  fit.  "Modes"  refers  to  the  number  of  modes  used  in  equation 
(28). 
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Figure  56:  See  Figure  55  for  a  description  of  the  models  being  compared. 


7  Summary  and  Conclusions 

One  long-standing  objection  to  the  Biot  model  is  that  its  predictions  for  com- 
pressional  wave  attenuation  are  not  in  accord  with  a  substantial  body  of  exper¬ 
imental  evidence  that  attenuation  in  unconsolidated  sediments  grows  linearly 
with  frequency.  The  measurements  made  in  the  SAX99  exercise  confirmed  this 
approximately  linear  relation,  at  least  at  frequencies  above  10  kHz.  The  use  of 
frequency-dependent  moduli  suggested  by  Chotiros  and  Isakson  [17]  (but  also  by 
Biot  [8])  in  lieu  of  the  frequency- independent  moduli  traditionally  used,  seems 
to  adequately  address  this  problem. 

There  are  two  reasons  to  challenge  the  Biot  model  on  theoretical  grounds. 
One  is  the  objection  of  Hickey  and  Sabatier  [22]  that  the  Biot  model  is  not  ap¬ 
plicable  to  heterogeneous  sediments  (Section  5.2.1).  This  criticism  is  valid,  but 
it  is  unknown  what  effect  this  deficiency  has  on  the  predictions  of  the  model. 
The  other  criticism  of  the  model  is  the  assumption  that  an  unconsolidated  sed¬ 
iment  can  be  modelled  as  a  Hookean  elastic  solid  of  constant  porosity  which 
encompasses  all  of  the  sand  grains.  As  discussed  in  Sections  5.2.2  and  5.2.3  this 
provides  a  feasible  explanation  of  the  poor  agreement  of  the  Biot  model  with  low 
frequency  data  in  the  SAX99  exercise,  at  least  when  the  measured  parameter 
values  are  used.  It  should  be  noted  however  that  the  presence  of  air  bubbles  in 
the  pore  space  provides  an  alternative  explanation  of  these  results,  one  which 
admits  the  possibility  of  air  bubble  resonance  which  might  be  of  some  utility 
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Figure  57:  Comparison  of  the  wave  speeds  predicted  by  the  four-mode  bubble- 
resonant  BICSQS  model  with  those  calculated  from  the  Kramers-Kronig  rela¬ 
tions.  See  the  text  for  further  explanation. 
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in  explaining  the  results  of  the  cross-hole  tomography  experiment  of  Turgut 
and  Yamamoto  [43]  discussed  in  Section  6.  As  noted  above,  however,  all  of 
the  postulated  mechanisms  make  predictions  too  similar  to  distinguish  among 
them. 

The  major  concern  about  the  Biot  model  that  remains  is  its  complexity 
and  the  number  of  parameters  that  are  important,  but  difficult,  or  at  least 
expensive  to  measure.  The  BICSQS  model  of  Chotiros  and  Isakson  does  not 
(disregarding  the  uninfluential  introduction  of  the  Poisson  ratio)  increase  the 
number  of  parameters  in  the  model,  however  it  does  replace  four  parameters 
Re  K£ ,  Im  Ky ,  Re  /i* ,  Im  p* ,  none  of  which  have  much  influence  on  compressional 
wave  speed  and  attenuation,  with  three  parameters,  the  gap  modulus  Ky  and 
the  relaxtaion  frequencies  fk  and  /;,  all  of  which  strongly  influence  the  BICSQS 
model’s  predictions  in  certain  frequency  ranges.  We  currently  have  no  means  of 
determining  these  parameters  other  than  from  measured  values  of  wave  speed 
and  attenuation  in  the  appropriate  frequency  range.  Likewise  the  distributed 
pore  size  version  of  Biot’s  viscosity  correction  factor  due  to  Yamamoto  and 
Turgut  [45]  should  be  used  in  lieu  of  Biot’s  original,  since  there  is  every  reason 
to  believe  that  pore  size  distribution  is  not  uniform.  This  again  replaces  a 
relatively  uninfluential  parameter  in  the  standard  Biot  model,  pore  size,  with 
the  standard  deviation  parameter  cr^,  which  as  indicated  in  Figure  14,  can  have 
a  substantial  effect  on  the  Biot  model’s  predictions,  at  least  if  ss  1  or  2  as  is 
conjectured. 

In  summary,  with  the  enhancements  suggested  by  Chotiros  and  Isakson  and 
Yamamoto  and  Turgut,  the  Biot  model  is  now  robust  enough  to  explain  all 
extant  measurements  of  compressional  wave  speed  and  attenuation  using  plau¬ 
sible  physical  mechanisms.  Its  predictive  value  is  still  in  question  because  of  the 
difficulty  of  determining  the  parameters  upon  which  it  depends. 
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